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XLV. Angular Distribution in (d, p) and (d, n) Reactions 


By A. B. Bratra,* Kun Huane,t R. Husy and H. C. Newns 
Department of Theoretical Physics, University of Liverpool ¢ 


[Received February 12, 1952] 


ABSTRACT 


The angular distribution of the particles emitted in (d, p) or (d, n) 
reactions is calculated, assuming a stripping process, by means of the 
Born approximation. The result resembles in most respects that of more 
complicated calculations by 8. T. Butler, but there are certain differences 
which are discussed in relation to experiments which ‘have been, or may 
be done. The particles are emitted mainly at fairly small angles to the 
incident direction, and the precise shape of the distribution is simply 
related to changes of angular momentum and parity occurring. The 
principal unknowns in the formulae are matrix elements, which may be 
estimated for sufficiently high energies by the statistical method. 


§1. INTRODUCTION 


Many measurements of the angular distribution in (d,n”) and (d, p) 
reactions have been made recently using deuterons of energy in the region 
8-15 Mev (e.g. Holt and Young (1950), Burrows et al. (1950), Hughes 
(1951), El Bedewi eé al. (1951 a, b), Gove (1951), Allen et al. (1951)). 

The most pronounced feature in general is the existence of one or 
more peaks at fairly small angles to the incident direction. A simple 
theory to account for these results was developed by K. Huang and 
A. B. Bhatia, and briefly reported on by R. Huby (1950). This theory 
is based on the stripping mechanism as proposed by Serber (1947) to 
account for (d, p) and (d, ) results (including the angular distribution) 
at much higher energies (~200 mev deuterons). At the lower deuteron 
energy now under consideration a simple Born approximation is applied 
to the stripping model, and a generally satisfactory form of angular 
distribution results. For outmoving neutrons or protons of definite 
energy, the form of the angular distribution is related in a simple way to 
the change of angular momentum and parity produced in the target 
nucleus by the capture of a nucleon from the incident deuteron. 

The use of the Born approximation for a nuclear process at the energies 
concerned is difficult to justify theoretically ; though as the form of the 
angular distribution turns out to be determined largely by such general 
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features as changes of angular momentum and parity, it may be hoped that 
the results are largely independent of the details of treatment, rather 
in the same way as the dispersion formula for nuclear reactions can be 
obtained by a simple second-order perturbation treatment (e.g. Bethe 
1937) as well as by more correct methods. 

Butler (1950, 1951) has also given a theory for the angular distribution 
in (d, p) and (d,n) reactions, based on the stripping mechanism. The 
mathematical technique used is not a perturbation one, but of a more 
penetrating kind involving the smooth joining up of wave functions for 
the captured nucleon at a radius 7, and thus allowing rigorously for the 
strong interaction of this nucleon with the nucleus. However, the theory 
appears still to rely on some approximations whose effect is difficult to 
assess, notably the assumption that there is a region near the nucleus 
(radius~ry) where one may neglect the interactions of the captured 
nucleon both with its former partner in the deuteron and with the target 
nucleus. It is believed that the development given in the general section 
(§ 3) of Butler’s (1951) paper contains some errors, though his final formula 
(34) is correct. For these reasons, and because of the complexity of the 
calculations of Butler, it is thought that the much simpler calculations 
of the Born approximation may be of interest, demonstrating as they 
do very clearly what actually happens in the reaction. 

The formula obtained with the Born approximation (18 below) is 
very similar to Butler’s in its general features. It contains a radial 
parameter R which is more arbitrary to fit correctly than the corresponding 
one (79) in Butler’s formula, but otherwise it usually gives equally good 
agreement with experiment. 

The angular distributions obtained with the two formulae differ most 
when the final nucleus is left in a state only just stable, or slightly unstable, 
against nucleon decay. The preliminary indications are that in such a 
case the Born formula may give the better agreement with experiment, 
but further experimental evidence is needed (see § 6). It should be noted 
that neither formula takes account of the Coulomb repulsion between 
the deuteron and the nucleus, so that it is necessary that the deuteron 
energy be sufficiently above the Coulomb barrier. Coulomb effects will 
presumably persist up to higher energies in the (d, ) than in the (d, :) 
reaction, because in the former the proton bearing the charge has to 
penetrate right to the nucleus, whereas in the latter the reaction can take 
place without the proton’s ever coming so close (cf. Oppenheimer and 
Phillips (1935)). 

The formula derived below (18) has been applied by H. C. Newns to 
experiments by El Bedewi et al. (1951 a, b) and El Bedewi (1952). 

The apparent success of the Born approximation for calculating angular 
distributions in the stripping reaction suggests its application also to other 
reactions which may depend on a related mechanism. This has already 
been tried for the inelastic scattering of deuterons, and results obtained 
consistent with experiment (Huby and Newns 1951). 
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The formula of Huang and Bhatia (18) is derived in § 2, discussed in § 3 
and compared with that of Butlerin§ 4. In § 5 it isshown how the formula 
can be used to obtain information about the state of the final nucleus. An 
example is discussed in § 6. The Appendix deals with results obtainable by 
recourse to the statistical theory of nuclei. 


§2. DERIVATION OF FoRMULA FOR THE CROSS-SECTION 


The discussions following will refer specifically to the (d, p) reaction, but 
the (d, ) reaction is embraced by interchanging the roles of neutron and 
proton throughout. 

In the stripping theory of the (d, ») reaction, it is assumed that the proton 
in the deuteron never comes into contact with the target nucleus: the 
neutron in the deuteron strikes the nucleus and is absorbed, while the 
proton continues its initial motion undisturbed, and emerges. We there- 
fore only take into account an interaction of the neutron with the nucleus. 
The Born approximation for the differential cross-section is then given by : 


c= age ies eee eee (1) 
where 
Taree, Op £) xp (G,)0 "PV (,, o,,€) x; (£) Ya (Fy —" ps Ins Op) 
eltka: ( atepi GrPdn donde, dC. BS 03h <== (2) 
The symbols have the meaning : 
Pee , MM, 
M,,=reduced mass of proton= M,+M,’ 
M*—reduced mass of deuteron= gu, 
Q M+; ; 


where M,, M,, Mz, M;, My, are the masses of proton, neutron, deuteron, 
initial nucleus and final nucleus respectively; k,=(M,Vv,)/h, 
k ,=(M7v,)/f, where v, is the velocity of the outgoing proton relative 
to the centre of gravity of the final nucleus, and v, is the velocity of the 
incoming deuteron in the laboratory system ; r,, r,, are the coordinates of 
the neutron and proton relative to the centre of gravity of the initial 
nucleus ; r, are the coordinates of the proton relative to the centre of 
gravity of the final nucleus (r,’=r,—(M,/M,)r,) ; o,,%,) are the spin 
coordinates of neutron and proton; & represents all the internal 
coordinates of the initial nucleus ; y,,, xg Xi X- are internal wave functions 
for the proton, deuteron, initial and final nuclei respectively ; V(r,,, o,, €) 
is the interaction operator between the neutron and the particles in the 
initial nucleus. 

We assume that the process of capture of the neutron takes place close 
to the surface of the initial nucleus, or, more precisely, in some narrow 
shell just outside the surface, where we can roughly imagine the attraction 
of the nucleus for the neutron to outweigh that of the proton in the 


2L2 
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deuteron. Let the mean radius of this shell be R, which ought to be 
approximately the radius of the nucleus plus the range of nuclear forces. 
Mathematically, this assumption implies that the major contribution 
to the integral (2) comes from a small range of 7, close to R, and so in 
the factors of the integrand which represent the incident wave, we may 
replace r=(r,,, 6,, ¢,) by R=(L, 6,,, ¢,) and obtain 


I=Jy,(o,)e ° ™» x al R— i Ge a, etka ok R+r,,)I/2 


{hyp (Ops On» E)V( Fao On> E)X iE)” dr, dé} dr, do, do, do»,,. Of ae (3) 


The expression in the brackets is a function of the angular coordinates 
and the spin of the neutron. It may thus be expanded as follows : 


Sxg (ns On» Vv Kno In» E\x(E)r? Ar, dé ; 
= a (f|V [i bn» Mn» bn Yy,™"(Ons Dn) Xin (Fn)> Se kay (4) 


Ins Mp» My 


where Y "is a spherical harmonic, and Xun (Fn) is a spin eigenfunction. 
The expansion coefficients are given by : 


CFV fe, ln» Mn» eo) ib ae on» é) Vics On> E)x;() Y,"" 


(9ns Pr) Xun (On) Wyn Ion dg. . - ~ (5) 
Substitution of (4) in (3) yields 
[= ect els Ins Mns Bn Klns Mins Bn PIA), + + + (8) 
where 
ny Mrs Mans P|EY=JVy,""" (Ons bn) Xin (Fn)Xp ple MPP’ 
XxXa(R—Kp, op, o,)elKa-(Rt tp? dr do, do,dw,. . (7) 


Now the components of angular momentum of the various internal wave 
functions Viz. uy, (Wg, 4; My are experimentally unspecified, and so the 
quantity ||? required in (1) is rather an average over the initial and sum 
over the final components, i.e., with the aid of (6) we modify |J |? to : 


1 
[Z Px 3(25. ~ 2 (Hy |[V [Mis Ens Mns Mn Uns Mins ns Mp [Ha ef 
3( Dit V una In» 


Ms ly 

: (8) 

Here the insertion of symbols such as y, in the Dirac brackets <|) 

indicates that internal wave functions Xu ete. with the specified angular 
momentum are to be used in the integrals (5) and (7). 

(8) can be simplified if the d-state component of the deuteron’s wave 


function is neglected, i.e. if yz is a product of a triplet spin factor and an 
s-state orbital factor :* 


Xug\ R— Le) On» Oy) =h( | R— Yy [Sy q(Ons Op): 


yy Ni 


* If the d admixture were allowed for, the result would probably be more com- 
plicated, and the selection rules governing J, (§3) would probably be less 
stringent. These effects are, however, expected to be of small magnitude. 
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hn» Mn Ens Mp LHad=9,"™{ bens Lp lea » pe CO'nd 13 > (9) 
where 
91,7 = JY." (Oy, bye “Fr'f( [R—r, [Jolt + dar,dey, . (10) 
and 
(Bp Ln |Ha = Tae pea (a re) Sugl@ n? Op) doy, doy. =aiice (11) 
After substituting (9) in (8), we can easily perform the summation over 
/, and wg with the result : 
i 
IP=>>—— 2 25 V peas bas Man Le lies 
P= 2(2), se Din Jon eae «Hs | le: 2» Mn Yin 
The summation over p,,, ~; and yz, can now be carried out with the aid of 
the group-theoretical properties of the wave functions concerned in the 
definition of the Dirac bracket expressions. The result is 


(12) 


i SN es hae lg, mn? (13) 
| lemeD ater Ree (OesAy See Tt) 
where 
A | (Hs lV [ois bn ’ Mn’ Lire) |e (14) 
mim 


91," can readily be evaluated from (10). If the principal axis for 
determining /,, is taken in the direction of 


M. 
=k,— >; iMaeetiz opel 
k=k, I, Na ene (15) 
1/2 
(note b=| (ba~ 5 gh +457 lglg sin? A 240 | , where @ is the angle 
between k,, and k,), then we Fe ; 
1/2 
Ii" =m, > 427( ly} 1)] 26 | k,, —$kcg |) (stn) Tie SiH), rear (16) 
where 
4a ats 
O(K)= = | init Speer.) (17) 
KJ} o 
Substitution of (16) into (13) and then into (1) yields 
ee ky a \ 1/2 2 eat 
A cane for Gis alae if we use oie the deuteron spatial wave function 
the approximate form : 
g \/2e-%r 19 
Hn= (5) 5, hee pee pene er) 
where «=1/h x (M,, x deuteron binding energy)*/ 
20-2310 om tee, oe eee ee (20) 
Then 
2(2aa)1/2 2(2rra)1/? 


G(|k,—2ka|)= (21) 


a+ |k,—3kgP  o®+(k,—$hq)?+2k,kq sin® 40° 


Dp 


49) A. B. Bhatia, Kun Huang, R. Huby and H. C. Newns on the 


§3. INTERPRETATION OF RESULTS 


The cross-section o (18) contains a sum of terms over a parameter 
|, which may be interpreted, in view of definitions (5) and (7), as the 
orbital angular momentum carried by the absorbed neutron into the 
nucleus. There are three important factors in the formula (18) : 


(i) Deuteron Factor G? (17 and 21) 

The proton, which initially has the mean momentum /k,/2 in the 
centre of mass frame, emerges with the momentum ik, The difference 
h(k,, —k,/2) represents the momentum of the proton in its motion inside 
the deuteron, at the instant of its release due to the removal of the 
neutron ; and the factor G?(|k,,—k,/2|) is proportional to the probability 
of this internal deuteron momentum. Equation (21) shows that this 
factor, as function of the angle 6 between k, and k,, has a forward 
maximum. The reason is that the larger the angle the larger the 
momentum difference which must be provided by the internal motion 
of the deuteron, and hence the smaller the probability. The factor 
G2 applies equally to all transitions with different J/,,. 

(ii) Neutron Centrifugal Factor [(7/2kR)¥2J iy ty( AR) 

From the last process, the neutron can be considered as issuing virtually 
from the deuteron with momentum #(k,—k,,) in the centre of mass frame. 
After its capture by the target nucleus, the neutron has the mean 
momentum —(M,/M,)hk,. The mean momentum interchanged between 
neutron and initial nucleus is thus hk=/[k,—(M,/M,)k,]. The factor 
[(n/2kR)*?J, .,(kR)]? is proportional to the mean probability that the 
neutron, travelling with this momentum in the frame in which the centre 
of gravity of the final nucleus is at rest, be found at the surface R of the 
nucleus, in a state in which the neutron’s total orbital angular momentum 
is U,. 

With the aid of (15), we see that the present spherical Bessel function 
factor is an oscillatory function of the angle 0, the oscillations decreasing 
as # increases. Provided |k;,—(M,/M,)k,|R is not too large, the factor 
has a maximum at the origin for J,=0, and for all other J,, a minimum 
at the origin, followed by a primary maximum which decreases in 
magnitude and moves to larger angles as J,, increases. The reason is 
that the larger the angular momentum J,, the greater is the momentum 
hk (and hence the angle @) which is necessary for the neutron to penetrate 
to the radius R. 


(iii) Nuclear Matrix Element A ;, (5 and 14) 


This is proportional to the probability of capture of a neutron which is 
at the surface of the nucleus with initial orbital angular momentum J,. We 
may note that A, is determined only by the properties of the initial 
and final nuclear levels, and by J,. When j,; and j, and the initial and 
final parities are given, (14) shows clearly that there will be selection 
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rules for the values of J, which yield non-zero values of A i, Lt follows 


from conventional group theoretical arguments that the selection rules 
are : 


(a) The final angular momentum j, must be obtainable by vector 
addition of the initial j,, the orbital angular momentum J, of the absorbed 
neutron, and the latter’s spin angular momentum 1/2. 

(b) J, must be even or odd, according as the parities of the initial and 
final levels are the same or different. These selection rules severely 
limit the number of terms appearing in (18), and can often bring it down 
to only one. 

It is not easy to make predictions about the magnitudes, absolute or 
relative, of the A, ’s. Generally, it will be expected that different final 
states of similar constitution will have A,,’s of roughly the same order 


of magnitude. For example, for any final state in which the captured 
neutron is approximately moving in an orbit outside the initial nucleus 
in its original state, there will be only one value of J, contributing 
appreciably, equal to the neutron’s orbital angular momentum, even 
though more J,’s might be permitted by the selection rules (a) and (6). 
The A,,’s for such final states might be considerably larger than those 
for most others ; and the relative values of the 4, ’s for different final 
states of this type would be largely governed by a statistical factor 
(2,+1). 

A rough estimate of the 4; ’s may be made by the statistical method, 
as used by Bethe (1938) for calculations of the cross-section integrated 
over all angles. It is shown in the appendix that, for highly excited 
final states unstable against neutron emission, this yields 

Wk, BR? (25,+1)e 
oe 40M pi (By) . Dt Sy OE cae (22) 
Here k, is the wave number of the neutron emitted if the final nucleus 
decays to the initial, i.e. 


2 ; 
< = —(Q of (d, p) reaction leading to final state f+ binding . (23) 
energy of deuteron). 


p;p(H,) is the density of levels of the final nucleus in the neighbourhood 
of the actual final level; counting only levels of the same angular 
momentum j, and parity as the level under consideration. . 

uw is 1 or 2, according as only one, or both, of the vectors of lengths 
(j,+4) and (j;—4) can be added vectorially to J,, to give a vector of 
length j, (otherwise, of course, 1.—0). 

The derivation of (22) only applies for k,,R>l,, but we may perhaps 
cautiously extrapolate it to smaller k&,, and even to imaginary 
k, (corresponding to a state of the final nucleus stable against neutron 
emission), if we replace a factor k,R in the numerator by unity—for an 


order of magnitude result. 
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It is also shown in the appendix how the use of (22) in the cross-section 
(18), for very high deuteron energy, reproduces essentially some of the 
results on angular distribution and energy spectrum obtained by Serber 
(1947), who did not analyse in detail the changes occurring in the target 
nucleus. For example, summation of the differential cross-section over 
many neighbouring final levels yields a result in which the combined 
effect of the centrifugal factors (spherical Bessel functions) is independent 
of angle, and the residual angular distribution is due to the deuteron 
factor alone. 


§ 4, COMPARISON WITH THE FoRMULA OF BUTLER (1951) 


For greatest ease of comparison of formula (18) with the formula (34) of 
Butler (1951), the latter may be rearranged somewhat, and the notation 
modified, to yield : 

iM Mik, 
ee eee (ae meal 
3 S(O aL) Ae (]k,—3k,]) 


1/2 
x 2; Ni eknto | Kuro, alto (sez) . Ji, 44(k1 0) 


In 


aler,\ 1/2 2 
+ Kj, 43(kn7o) ( 5} ) J,,-a(kro) | . . : ° . é . . (24) 


Here @’ is a function (deuteron factor) similar to G in (21), but differing 
a little because Butler used a more accurate form of wave function for 
the deuteron ground state than (19). 7)isanuclear radius. «, is defined by 


Kyi — tk; ee ere 


(cf. eqn. (23)), and is thus real if the final nucleus is stable against neutron 
decay. The case that the final nucleus is formed in a state unstable 
against neutron decay, and hence « is imaginary, does not have to be treated 
separately in the Born calculations but it must be in Butler’s. While 
Butler does not discuss this case fully, it is understood that the cross- 
section is still given by (24) provided the actual imaginary value of 
x, 18 inserted, and the square bracket is replaced by a modulus. 

Remembering that the only variables containing the angle in (24) are 
k (eqn. 15) and {k,—3k,|, we see that the only essential difference in 
the angular distribution between (24) and (18) lies in adding to the spherical 
Bessel function (7/2kr9)¥2J i,+4(479) another term proportional to 
(mkro/2)'J, (kro). The effect of this is usually rather small, and either 
formula can usually be made to fit an experimental curve equally well, if 
suitable values of 7) and R are chosen, the latter always having to be rather 
larger than the former. : 

The difference depends mainly on «,. It is very small for |x, | large 
(state of final nucleus very stable or very unstable), and greatest for 
|x, | small (state of final nucleus only just stable or just unstable), 
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because in this case the second term in the square bracket in (24) takes 
control. The effects ensuing will be illustrated in an example in § 6. 
The conclusion is that Butler’s formula behaves in a ‘ special ’ or ‘ singular ‘ 
manner when x,—0, while the Born formula does not, and the precise 
behaviour of the former will depend rather sensitively on the value of 
x, When small. It does not seem clear whether this singular behaviour 
is to be expected physically, and it would thus be of interest to make 
careful comparisons with experiment in this region. 


§5. APPLICATION TO DETERMINATION OF STATES 


It is to be expected that the formula (18) will only agree with 
experimental curves in the region of small angles where the cross-section 
is large ; since the constructive interference of many partial waves of 
definite orbital angular momentum for the proton and deuteron, which 
is really responsible for the maximum cross-section, should be more 
reliable than the details of the destructive interference which makes 
the cross-section small at large angles; and, in addition, there will 
probably be a background of neutrons produced by the competing 
process involving formation and decay of a compound nucleus, and 
these will tend to mask the stripping cross-section where it is small.* 

To compare (18) with an experimental curve, it is necessary to know 
which A; ’s are present, and their value. For certain initial and final 


states more than one J, may be permitted by the selection rules of 
§ 3 (iii)—though the values permitted must be either all even or all odd. 
However, in practice, it has been found possible in most cases to obtain 
a satisfactory fit with one /, only. When this has been found, if the parity 
and spin of the initial state are known then the selection rules fix the 
parity of the final level uniquely, and its spin j, to within the limits 
permitted by vector addition of j,; plus /, plus 4. 

To find J,, the experimental curve should first be divided by the 
theoretical deuteron factor G? (see 21), so that the result is to be 
proportional simply to the sum in (18). A first indication of J,, can be 
obtained by inspection of this curve at small angles : a forward maximum 
will normally indicate that J,—0 is present, and a forward minimum 
that it is not (though these rules are not absolute). We can calculate the 
several spherical Bessel function terms in (18) and find which gives the 


* Cohen and Falk (1951) put forward experimental evidence from the energy 
spectrum in (d, n) reactions that the forward scattering is mainly due to stripping 
and that at 90° to compound nucleus. It is usually assumed that the angular 
distribution of particles from the compound nucleus will be fairly isotropic, 
but Wolfenstein’s theory (1951) predicts peaking in the forward and backward. 
directions, with symmetry about 90°. Indeed, it seems that Wolfenstein’s 
theory could lead in extreme cases to a forward peak of sharpness comparable to 
that due on the stripping theory to the ‘ centrifugal ’ factor § 3 (ii). However, 
it should be possible to tell whether such a forward peak is due to the stripping 
or the compound nucleus process, because the latter should be accompanied by 
a backward peak, and the former not. 
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best fit with experiment, if we take a suitable value of k—and this last 
point, to choose the right R, is a major consideration. Butler (1951) 
found that with his formula (24) he usually obtained good and consistent 
agreement using for his radius 7) an accepted value for the radius of the 
nucleus plus range of nuclear forces e.g. as given in Gamow and 
Critchfield (1949, p. 11): 

radius =(1-7-+1-22x Al/)x10-%cem; . . . . (26) 
and on the average for light nuclei he took 5-0x10~% cm. However, 
to obtain good and consistent results with the Born formula, a larger 
value of R is found to be needed, on the average 2-0 x 10~1!8 cm greater 
than the Gamow value (26). 

The procedure is thus to find what value of /,,, together with a value of 
R not very different from that just mentioned, will give the best 
agreement with the experimental curve. The answer has in practice 
almost always turned out to be unique, for a given curve. 

The large value of Rk required is hard to justify physically, and is in 
contrast to the reasonable value of the 7, required in Butler’s formula. 
Some investigation of the role of & or ry has therefore been made. 

One hopes that the dependence of the angular distribution on & or 
ro will be insensitive, so that their exact choice is not very important. 
To obtain a measure of the sensitivity, we follow up a suggestion of 
Professor R. E. Peierls. Take the cross-section for /,—=0 according to 
(18) and let ky be the value of & at which its first zero occurs. Then we 
may use as the gauge of sensitivity the quantity (R/k,)(dk)/dR). Since, 
in fact, kyocR, this is just 1. The corresponding quantity (79/k)(dko/dr9) 
for Butler’s formula (24) les always between 0-78 and 1. This shows 
that the sensitivities of the two formulae are nearly the same. 

Next we may analyse the statistics of a number of reactions for which 
both (18) and (24) have been fitted to the experimental curves. In each 
case the best values of R and 7, to give a fit have been found, and compared 
to the Gamow radius (26). The results for 10 curves are that : 


for Born formula (18), R=Gamow radius+2-0 x 10-3-+-0-50, 
for Butler’s formula (24), 7>=Gamow radius+0-36 x 10-3-.0-44. |” (2 


(The limits are * probable errors “0-67 x R.M.S. deviation from mean.) 
Thus the spread of the best values of the radius is about the same for both 
the Born and Butler’s theories ; but on the average Butler’s radius comes 
out closer to the physically reasonable Gamow radius. Both methods 
appear in practice to be about equally reliable for finding the 
characteristics of the final state. 


§ 6. EXAMPLE 
We shall discuss only the results of the reactions O18(d, n)F!" and 
OM(d, n)F™*, which have been obtained by El Bedewi et al. (1951 a) 
using 8 Mey deuterons. The Q’s of the respective reactions are —1-36 Mev 
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and —1-89 Mey, yielding, in accordance with (25) and (23) the respective 
values for x,t of 0-196 x 1033 em-! and 0-119 1018 em-!, the latter in 
particular being unusually small. 

Considering first the reaction forming F1" in the ground state, we find 
that both eqns. (18) and (24) require Lt =2 for a good fit. The radius 
R giving the best fit by Born approximation is 6-25 10-13 em (note 
that the Gamow radius from (26) is 4-8 10-18 em, and compare with 
(27)) ; while with Butler’s formula, a good result is obtained with the 
radius 5-0 10-18 em which he recommends on the average. The results 
are shown in fig. 1. The insensitivity of the results to the value of the 


Fig. 1 


DIFFERENTIAL CROSS-SECTION 
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Angle (Centre of Mass System). 
Reaction O16 (d,n) F!” with 8mev deuterons. Q=—1-36 Mev. 
Points are measured by El Bedewi et al., 1951 a. 
Full curve : theory according to Born approxn., with Rk =6-25 x 10~-¥ em. 
Broken curve: theory according to Butler, with ry=5 x10-™ cm. 


radius is shown by the curves of fig. 2, which are drawn for radii r, and 
R each increased by 110-1 cm from the values mentioned above : 
the fit is still quite good. The reliability of the value of J, is indicated 


by the fact that, taking the Born approximation, for instance, to obtain 
eengepeneey teedet a So te ES ei 

+ Note that, as this is a (d, ) reaction, the roles of neutron and proton are 
interchanged from those in the preceding discussions. 


496 A. B. Bhatia, Kun Huang, R. Huby and H. C. Newns on the 


a good fit assuming /,=1 or 3, we should require respective values of 
R of 3-9X10-23 and 8410-3 em, neither of which conforms at all 
well to (27). 

The conclusion to be drawn from the selection rules of § 3 (iii) is that, 
since the ground state of 016 has j,;=0, even parity, and since /,=2, then | 
the ground state of F!” has the spin j,=3/2 or 5/2 and even parity. 

For the first excited state, both formulae (18) and (24) require 
l,=0. The best values of & and ry for a fit are respectively 7:25 10~-# 
and 4-2 x 10-13 em, and the curves for these are shown in fig. 3. Increasing 


Fig. 2 


DIFFERENTIAL CROSS-SECTION 
(ARBITRARY UNITS) 
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Angle (Centre of Mass System). 


Reaction O18 (d,n) F1? with 8 mev deuterons. @—=—1-36 Mev. 
Points are measured by El Bedewi et al., 1951 a. 

Full curve : theory according to Born approxn., with R=7-25 x 10-18 em. 
Broken curve: theory according to Butler, with r7=6 x 10-13 em. 


ro to the ‘average’ value quoted by Butler, 5-0 10-™ em, gives poor 
agreement, as shown in fig. 4, whereas an increase of R in the formula 
(18) to 8:05 x 10~!8 does not cause so much deterioration. The reason 
is the ‘ singular’ behaviour of Butler’s formula for small « (see § 4), two 
of the chief features of which can be seen in this example, viz. ‘(i) a 
reduction of the value of kr) at which the first zero of the cross-section 
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occurs (tending to require unusually small values of r,) and (ii) unusually 
large secondary maxima at large angles.* 

The result 1,=0, together with the known spin and parity of O15, leads 
to the result j,=4, and even parity, for the first excited state of F217. 


Fig. 3 


DIFFERENTIAL CROSS-SECTION 
(ARBITRARY UNITS) 


en One 0s AQ S02. 60% 70? 
Angle (Centre of Mass System). 

Reaction O16 (d, n) F1’* with 8 Mev deuterons. Q@=—1-89 Mev. 

Points are measured by El Bedewi et al., 1951 a. 

Full curve : theory according to Born approxn., with & =7-25 x 10~!8 cm. 

Broken curve: theory according to Butler, with 7)=4-2 x 10-¥ cm. 


APPENDIX 


Use of the Statistical Theory of Nuclez 
(Following Bethe, 1938) 

We shall express the matrix element 4, (14) in terms of the partial 
width for decay of the final nucleus with neutron emission, leaving the 
initial nucleus in its initial state (assuming this to be energetically possible). 


* Although there is as yet not much evidence under these ‘singular’ 
conditions, it is of some interest that in one other case, CO! (d, n) N13*, (see El 
Bedewi et al. 1951 b) forming the final nucleus in its first excited state, which is 
slightly unstable against neutron decay («,=—tx0-14x 10% cm™'), it does 
not prove possible to obtain a really good fit for Butler’s formula with any 
value of ry: the nearest fit is got with 7) =3-7 x 10-18 cm, whereas 77 =5 x 10-1? em 


is very poor. 
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The partial width can then in turn be estimated in terms of a sticking 


probability. 
The formula for the partial width for the decay process in which the 


neutron is emitted with angular momentum 1, is, by perturbation 
theory : 


4h, 2 be * My * 4 ys ¥ 
UD ryote h?(2j,-+1) Migs Mf I Xj E) Xu (Fn) Tp us (6, by)(7/2k,7 nd Si 4 (Kn! a) 


Lyn, My, 


x Vras Ons Xu as Ops ddr: dog dE Aa eee ee) 


Fig. 4 


DIFFERENTIAL CROSS-SECTION 
(ARBITRARY UNITS 


fe) 109299 20° S020 AO Le SO Ome GO on O= 
Angle (Centre of Mass System). 
Reaction O18 (d, n) F17* with 8 Mev. deuterons. Q—=—1-89 mey. 
Points are measured by El Bedewi et al., 1951 a. 
Full curve : theory according to Born approxn., with R=8-05 x 10-13 em. 
Broken curve: theory according to Butler, with ry=5 x 10-13 cm. 


(The notations are those of § 2 ; k,, was defined in (23).) If we make the 
assumption, as in § 2, that the interaction is effective only at the radius 
fk, so that r, in the spherical Bessel function (28) may be replaced by 
R, then with the help of (5) and (14) we obtain at once : 


4k, A, a \1/2 2 
Dis FED | SER) Jip sally) | sfpesek aay 
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If k,R is appreciably larger than (J,,-+4), the asymptotic form of the 
Bessel anction may be used, yielding on an average over a range of k, 


Mai, 


Liens = h?(2),+ 1k, R (30) 


Now by the principle of detailed balance the partial width is expressible 
in terms of the capture probability for a neutron of wave number 
colliding with the initial nucleus (see e.g. L. Wolfenstein, 1951) 


| An AV GAD EF, s sity 31) 
delas dp (297+ 1) p; (Hs) : coer Geshe! Rete 


p;,(H;) was defined below (23). Figs ipl, IS the sticking probability 
for the neutron, if its orbital angular momentum is 1,, to combine with 
the initial nucleus j; to form a compound, nucleus with the angular 
momentum j, and proper parity.* 

Now if we assume that the sticking probability does not depend 
explicitly on the final angular momentum j,, but is determined only by 
statistical factors, we may write 


Ky 


ens (2+ 1)u Es, 
Ifijisln 2(21, £1)(27,41)’ 


where » was defined below (23), and ¢,.;, is the sticking probability 


(32) 


irrespective of final state, which to the present accuracy may be taken 
as 1. Substituting (32) in (31) and thence in (30), we obtain the final 
formula (22) for A, . 

We now investigate the cross-section obtained in measuring many 
proton groups together. We find first the differential cross-section 
do(@) for all protons in the energy range dH, scattered at an angle 
8, from (18) : 


MMi k cra 2 
do(8)= Inj Eee he a lice —3k,|). dz, feos (sen) Fiyss(tR) | 
x pj (Es) 
(substituting from (22) and summing over j;) 
G?(|k,—4k,q|)R? Mykyk, OTR 1/2 2 
ee oe eed 1 Pde po (2h) ER Ti iy(kR) |. 
BG Hg |) RF Malek Vay au, COPEL ue nee 
—_ ee be — dH A ee (33) 


since the sum over J, is just unity. (In the above we have put M,*=M,, 
etc.) The only angular dependence here is through the deuteron factor 


* The sticking probability as here used has to contain also any inhibiting 
effects of barrier penetration. 
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GS (21) shows that, for high deuteron energy, this has a sharp peak, as 
function of k,, when k,=}k,. For important cases, then, we may put 
in (33) 


ka : 
k,=k,= coe and obtain 


RM, kg 
Sa7h? 
We can integrate (34) over all proton energies H,,, with the aid of (21) ; to 
obtain a total angular distribution o;,,,(@) (after approximating for 

large ky) : 


do(0)=G?(|k,—$ka|) di Sa eee sa 


RAke Ao. 
wt) = CETTE E sin? OPE [oy ei eee eeh 


The energy spectrum do for all protons in the range dH, is obtained 
similarly by integrating (33) over all angles 0 with the aid of (21) : 
ee 16R2M jak, k,dk , y 
hk ala? +(ky—gkq)*(ky+ oka)? 
The results (35) and (36) agree essentially with the eqns. (11) and (15) of 
Serber (1947) for the angular distribution and spectrum assuming a 
‘transparent’ nucleus, i.e. neglecting entirely the interaction of the 
proton with the target nucleus (only in Serber’s paper it would be the 
neutron, not proton, since he discusses (d, ”) reactions). 
The total cross-section o;,, for all angles and all proton energies is 
obtained by integrating either (35) over angles or (36) over #, : 


(36) 


Opp 27’, ° . . . . . . (37) 
which is of the right order of magnitude, though too large. 
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SUMMARY 


In this paper the theory of the stability of a layer of fluid heated below 
which has been developed by Rayleigh, Jeffreys, Pellew and Southwell 
is extended to the case when the fluid considered is an electrical conductor 
and an external magnetic field is impressed on the fluid. The problem 
is first investigated on the assumption that the principle of the exchange 
of stabilities (i.e. the equations governing marginal stability are those 
which are obtained from the general time dependent equations by setting 
o/dt=0) is valid. A differential equation of order six for the normal 
component of the velocity is derived. Suitable boundary conditions are 
formulated ; they depend on whether the bounding surfaces are free or 
rigid. The case when the magnetic field and gravity act in the same 
direction is investigated in some detail. A variational procedure for 
solving the relevant equations satisfying the necessary boundary conditions 
and determining the critical Rayleigh numbers for the onset of convection 
are described. Tables of the critical Rayleigh numbers are provided for 
the three cases (i) both bounding surfaces free, (ii) both bounding surfaces 
rigid and (iii) one bounding surface free and the other rigid. Finally, 
the principle of the exchange of stabilities is examined and it is shown 
that it is valid if 4zuox <1, where p is the magnetic permeability, o the 
electrical conductivity and « the thermometric conductivity. This 
condition is satisfied under most terrestial conditions; but under 
astrophysical conditions it is not. In the latter case, it is shown that 
instability can arise either as cellular convection or as over-stability 
(i.e. by oscillations of increasing amplitude) depending on the magnitude 
of the velocity of the magneto-hydrodynamic wave (=(uH?/4zp,)"”). 


§1. INTRODUCTION 


THE manner of the onset of convection in an incompressible fluid by 
thermal instability has been the subject of investigations by Rayleigh 
(1916), Jeffreys (1926, 1928), Low (1929), Pellew and Southwell (1940) 
and others. The principal results established by these investigations 
are (i) that a layer of fluid heated below first becomes unstable when the 
‘Rayleigh number ; 
pals eae cee, 


K 


*Communicated by the Author. 
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(where g denotes the acceleration due to gravity, d the depth of the layer, 
B=—|B| is the adverse temperature gradient which is maintained and 
a, « and v are the coefficients of volume expansion, thermometric 
conductivity and kinematic viscosity, respectively) exceeds a certain, 
determinate, critical value of the order of a thousand and (ii) that the 
motions which ensue on surpassing the critical value must have a cellular 
pattern in agreement with what had been demonstrated, experimentally, 
by Bénard (1900, 1901). The exact value of the Rayleigh number at 
which instability first sets in depends on boundary conditions such as 
whether the plane surfaces between which the layer of fluid is confined 
are free or rigid; and the values which have been derived theoretically 
have been confirmed by the experiments of Schmidt and Milverton 
(1935) (in this connection see also Sutton 1950). 

While the criterion for instability established by the methods of 
Rayleigh and Jeffreys is derived explicitly for an incompressible fluid, 
Jeffreys (1930) has shown that the same criterion can also be applied to 
compressible fluids provided the density does not vary greatly within 
the system and provided also we interpret | 8 | as the temperature gradient 
in excess of the adiabatic. Moreover, it follows from the discussion of 
Rayleigh and Jeffreys that when treating an incompressible fluid in which 
the temperature is variable, we can often neglect the variations of density 
(caused by thermal expansion) except in so far as they modify gravity ; 
this is a method of approximation which apparently goes back to 
Boussinesq (1903). 

Now the occurrence of magnetic fields in a number of geophysical and 
astrophysical problems in which convection is presumed to play a role 
makes it of interest to examine the effect of an external magnetic field 
on the onset of convection by thermal instability. It would appear that 
a profitable beginning in this direction may be made by re-examining 
the Rayleigh—Jeffreys problem for the case when the fluid considered is 
an electrical conductor and an external magnetic field is impressed on the 
fluid. We shall then have to add Maxwell’s equations to the equations 
of hydrodynamics and heat conduction. However, it will appear that 
the addition of Maxwell’s equations to the other equations of the 
problem does not introduce any essential element of difficulty ; indeed, 
we shall see that the problem can be solved to the same degree of 
completeness that the classical Rayleigh—Jeffreys problem has been 
solved. 

After this paper was completed for the press, the December 1951 issue 
of the Philosophical Magazine was received at the Yerkes Observatory. 
This issue contains a paper by W. B. Thompson (1951) in which the 
problem to which the present paper is devoted is also considered. But a 
comparison of the two papers will show that except for the basic idea 
there is very little overlapping between them: thus in the present paper 
the mathematical problem is treated with due attention to boundary 
conditions (which are effectively ignored by Thompson) and tables of 
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the critical Rayleigh numbers for a wide range of possibilities are 
presented. The limits of validity of the principle of ‘exchange of 
stabilities ° are discussed in both papers; but it would appear that the 
present paper goes somewhat further than Thompson’s: in particular 
we include the effect of viscosity which is ignored in this connection by 
Thompson but which turns out to be rather important for astrophysical 
applications. 


§2. THE EQuaTIONS OF THE PROBLEM 


The equations of motion, of Maxwell, and of heat conduction appropriate 
to the problem on hand are: 


0 
Pa; +p(u. grad)u—pJ x H=—grad p+pvVu+p grad V,  . (2) 


curl H=47J, curl poe td ; 
eee! (3) 
div H=0, J=o(E+puxH) 
oT 
and a 7 (u OES iN LG Bee ol ye eee) 


where p denotes the density, p the pressure, 7’ the temperature, V the 
gravitational potential, E the intensity of the electric field, H the intensity 
of the magnetic field, J the current density, » the magnetic permeability 
and vy, « and o are the coefficients of kinematic viscosity, thermometric 
conductivity and electrical conductivity, respectively. It will be 
noticed that in writing Maxwell’s equations we have ignored the terms 
in the displacement currents; as is well known this is a permissible 
approximation in problems of magneto-hydrodynamics. (It may also 
be noted here that we are adopting the electromagnetic system of units.) 

Substituting curl H/47 for J in the equation of motion (2) and making 
use of the solenoidal property of H, we can reduce it to the form 


Ou; 7) pu 0 ete 0 e| HP 
pF toatl Fae Hilly — a (p a 
oV 
+ pvV2u +p On; 4 ed Wie Mae) EN «eae (5) 


where here and in the sequel summation over repeated indices is to be 
understood. Following Rayleigh we shall allow in eqns. (5) for the 
variation of density only in so far as it modifies the effect of the external 
field. Thus, in eqn. (5) we replace p which occurs in front of 0V/dx; by 


p=po(l—a«AT) ° ° . ° . . . . (6) 
where « denotes the coefficient of volume expansion, py the density 


corresponding to a certain mean temperature 7) and A7’ is the deviation 
of the local temperature from 7’): 


AT—T—T,; : ns ‘ 5 é ci 5 : i) 


2M2 
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and we regard p occurring elsewhere in eqn. (5) as a constant equal to py. 
On these assumptions eqn. (5) becomes 
0 OV 


Ou, — 0 fs ov 
sae! hs Ae (uuj)— Trp dai iA; =vV*u; “AT — e 


d ul HP 
—3,(% a ee -v). aS} 


With the variation of density due to thermal expansion allowed for in 
this manner, we treat u; from now on as a solenoidal vector : 


Ou; 
ea Mere, ars fy 
When the external field is that due to the action of gravity, we can write 
oV 
ae eae and V=—gA,z;,, (10) 


where A is a unit vector in the direction of the vertical. With the 
foregoing substitution, eqn. (8) becomes 


iad Aig arpa aa HNL EEE. ones TL SNE iy 
att da, 5 amp et es uryATr, 


a/p  ulHP 
~3,(4+ Baek “tb GAyty |G te es 


where we have written 


y=ga. ETP Re ar ho, ote: 
Turning to the equation of heat conduction (4), we shall suppose that 
a constant mean temperature gradient B=—|f| is maintained in the 
direction of A by an external agency and that we can write 
T=T,+BAjv;+0, AT=BA,x,+0 eM easy oy eS 
where 7’) is a constant (ef. eqn. (7)) and 
aq 
f= On; X;; . . ° . . . ° . (14) 


and @ is the deviation of the temperature from its local mean value 
T')+PA,x;. Inserting (13) in eqns. (4) and (11) and remembering that 
H, is solenoidal, we obtain 


Ou, 0 be 0H ; Ow 
ee Se Ty oe .—— —ypV2u, —— 
ot a Ox od J Arr py Ht 0x; vV uty, Ox; (19) 
00 06 
and ‘at TBAjus Uj Ox; =KV°, e . e e . ‘ (16) 
where, for brevity, we have written 
p eH 
oa eres + gAjt;—bByAjAa ee, . . . . 17) 
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Considering next Maxwell’s eqns. (3), we obtain for H the equation 


0H J 1 
ers =— curl (2 —ux) =p curl (ux H)— - curl J 
i) 
Curl (ucla VA Bae Ste! aretab dow. (18) 
dro 
Hence Be 2 
: or —curl (ux H)=7V?H, eens te Bee OR), 
1 
where = 
ere Career (20) 


Using the solenoidal property of u; and H,, we can rewrite eqn. (19) 
alternatively in the form 
0H, oH , Ou; 
—— ; t . Z — 2 ° 
aT siers Wi GN mnt ors te (ok) 
If a constant magnetic field H=|H| is maintained, externally, in the 
direction of a unit vector t, we may write 
H =H7-+h,, . . . . ° . ° . (22) 
where h,; is the perturbation in the magnetic field caused by the motions 
and the induced currents in the fluid. 
We shall now suppose that w,, @ and h, are small quantities of the 
first order and that we can ignore products and squares of them. On this 
approximation, eqns. (15), (16) and (21) reduce to 


du; pH dh; Ow 
Pee eae, 8) 
a0 
ee PERV eee Oe Ae (is. (2) 
a thee 

and am Pe ea 


Taking the divergence of eqn. (23) and remembering that uw, and h; 
are both solenoidal, we find that 
: 06 
pee ewes 
Rrra A re oti oa. eB, (26) 
In the further discussion of the foregoing equations, we shall find it 
convenient to work with the velocity, w, and the intensity of the magnetic 
field, 4, in the direction of A. Thus, letting 


w=u,A; and h=hA; ° ° ° . . . . (27) 
we readily derive from eqns. (23) to (25) the following set of equations : 
ow pH oh ; Cm 
Ot = Aor po Tj On, +y0-+vV w—A,; Ox,’ Pron ade eS (28) 

6 
x = —Bw+kv"d, ok 0 ge OD) 


oh ow 
and AF =AHr, an; +nV*h . ahs hal by te og 10 eR ERG ries sane (30) 
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These equations together with eqn. (26) represent the basic equations 
of this theory. 

As stated in § 1, we shall consider a horizontal layer of fluid of depth, 
d, confined between two parallel planes at z=0 and z=d. The conditions 
to be satisfied on these bounding surfaces will depend on whether they 
are free or rigid. But in all cases, we must require that 


6=w=0 for z=0 and z=d, ater Ie) 


since (by assumption) the temperatures on these surfaces are maintained 
constant and there can also be no motion normal to them. Additional 
boundary conditions follow from the equation of continuity (9) and 
depend on whether the bounding surfaces are free or rigid: they are 
(cf. Pellew and Southwell 1940) : 


02w 


re =()- ona free surface i. 70 aa ee eee 
dw . . 
and a =0 onarigid surface. 0 Sie A 


§3. THE EQUATIONS GOVERNING MARGINAL STABILITY 


If the principle of the exchange of stabilities as described by Jeffreys 
(1926) is valid also for the problem under consideration, then we can write 
down at once the equations governing the situation in marginal stability 
(i.e. one on the verge of stability) by setting 0/¢t=0 in eqns. (28) to (30). 
In § 8 we shall examine the validity of this principle for the problem on 
hand; and it will appear from the discussion there that the principle 
will be valid under most terrestial conditions in which we shall want to 
apply the theory. We shall therefore accept the principle for the present, 
postponing to § 8 an examination of the limits of its validity. 

Setting 0/dt=0 in eqns. (28) to (30) in accordance with the principle 
of the exchange of stabilities, we obtain the following equations governing 
marginal stability : 


do wi ah 


vV2w= A; Ox, ae ae Spy Tj x, > (34) 
xV20=Pw, dosnt Op coe. ee 

dw 

2h —— 
7V7h= AT 5,,? Sauer eerie Coons iit 

00 
27 — 

and V bale irr of eA Si ees Ae aes re! Te eee (37) 


Operating eqn. (34) by V2 and making use of eqns. (20), (36) and (37) 
we obtain 


L) 


ie 076 2772 aap 
viw= 2 (ra ie on dha bos . (38) 


iA; a a 
v\ "4 On, pov ' 1 Ox, 0x, 
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Operating this equation once again by V? and making use of eqn. (35), 
we have 
Ow tele 0 


enn WE Ae ee ee er 
WwW om V2w NUE -- oa "149 Ou, VO a9) 


It is also evident that 0 satisfies an equation of exactly the same form as w. 


$4. THe Repvuction To A BounpaRy VALUE PROBLEM WHEN 
BoTH H anD g AcT IN THE SAME DIRECTION 


In some ways the case of greatest interest is when the external magnetic 
field and gravity both act in the same direction. If this common direction 
be chosen as the z-direction, then 


heat OS 0c Ea eee eo 8h gt (40) 
and eqns. (38) and (39) of the preceding section reduce to 
ET i a 
ies “(aes oy? poy. O27 Pats veer: 
and 
jy LEA A or Wi a 
SAI = or (Fe aye ae ae ¥ ra sl ee (42) 


In analogy with the classical Rayleigh—Jeffreys problem, we may expect 
that in the present case also, convection when it first sets in will have a 
cellular pattern with the cell walls being surfaces of symmetry. We may 
then separate the variables x and y from z by assuming that w and 0 
have the forms 


w=f (a, y)W(z) and 0=f(x, y)@(z), . - . . (48) 
where 
Cains eal 
d? (5 + 53) = —a*f, Oe fC. = ON Oa (44) 


_and a? is a certain ‘ characteristic number ’ determined by the conditions 
on the cell walls. And if as we have assumed, the cell walls are surfaces 
of symmetry, we must require that the normal gradient dOf/én vanishes 
on the cell walls. As Pellew and Southwell (1940) have shown this 
boundary condition on f ensures the reality of a. 
The solution for f when the boundary of the cell is a rectangle of 
sides LZ, and L, is, clearly 
{—=Ay. eos “= 00s 7 seirelbeb at sha sie (ED) 
where m and n are arbitrary even integers and 4A,,,, is a constant; the 
corresponding relation between a, L,, L,, and the depth of the cell, d, is 


ihe ea! 
adn? (Ta + 7) . a z ‘ y 46) 
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The solution of eqn. (44) when the boundary of the cell is a regular 
hexagon with a side LZ has been discovered by Christopherson (1940). 
His solution is 
2n 2na 4niy 
Phy 00s = (o/3%+-Yy)+ COS SF (/3x—y) + Cos SF \ pa) paerleay 
where f, is a constant and n is an integer; and the equation relating 
a, L and d is (cf. Pellew and Southwell 1940) 


a= (Fr) Su ts eee aoe 


Returning to eqns. (41) and (42), substituting for w and @ in accordance 
with eqns. (43) and (44), and letting 


d 
C=7zd sand 1 de? (49) 
we have 
[(D?—a??$—QD"|W =" a0, oy o> Se 
and 
| (D?—a?)[(D?—a?)?—QD?|W =—a2RW, ere we iiss: 
where 
2772 
eee: and palBly ge, Saker So Se 
Pov KV 


The boundary conditions with respect to which eqns. (50) and (51) 
must be solved are (cf. eqns. (31) to (33)) : 
@=W=0. for C=0 and) 27. ee ees 
and 


DW=0 ona rigid surface, ] 
(54) 
D?>W=0 ona free surface. 


In view of eqn. (50), the boundary conditions (53) are equivalent to 
W=[(D?—a*}?}—QD*]W—0 for ¢=0 and) 3 3165) 


On a free surface the boundary conditions (54) and (55) can be combined 
to give 


W=D?W=D4W=0 (ona free surface). . . . (56) 


Equations (54) and (55) provide six boundary conditions and the 
requirement that a solution of eqn. (51) satisfies these conditions will 
determine for any given a? a sequence of possible values for R. For, 
eqn. (51) being a linear equation with constant coefficients of order six, 
its general solution must be of the form 


6 
w= PX A, e%®, . . . . . . . (57) 


t=1 


where the q;'s, occuring in pairs, are the roots of the characteristic 
equation 


(q?—a*)[(¢?—a?)?—Qq?]=—a2R, Roy as GI) 
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and the A,’s are constants of integration. The six boundary conditions 
(three at ¢=0 and three at ¢=1) will lead to a set of six linear 
homogeneous equations for the A ;’s and the condition that the determinant 
of this system vanishes (so that we may have a non-trivial solution) will 
determine a sequence of possible values for R. Among these possible 
values (for a given a?) there will be a lowest value; the minimum of 
these lowest values (as a function of a2) will finally specify the critical 
Rayleigh number, f#,, at which convection, through instability, will first 
set in. 

While there is, of course, no difficulty of principle in carrying out the 
procedure outlined in the preceding paragraph, it is evident that the 
actual carrying out of the procedure for various values of Q and boundary 
conditions will involve a prohibitive amount of work. Thus, even for 
the classical Rayleigh—Jeffreys problem the precise determination of R, 
for the two standard cases (i) both bounding surfaces rigid and (ii) one 
bounding surface rigid and the other free, was completed only in 1940 
by Pellew and Southwell. However, it will appear that the extension of 
the variational principle established by Pellew and Southwell for the 
Rayleigh—Jeffreys problem to the present one leads to a convenient 
practical procedure for determining R, without an undue amount of 
labour. The variational principle is developed in the following section. 


§5. A VARIATIONAL PRINCIPLE FOR DETERMINING fk, WHEN 
BoTH H AND g ACT IN THE SAME DIRECTION 


First, we shall obtain an integral expression for Rk. Letting 
F=(|(D?—a)?— QD? |W, AN thay eh age a (9) 


we shall rewrite the differential equation governing W in the form 
(cf. eqn. (51)) 
(D?—a?)F=—@RW. Sia Nbewie ei? oa regs (OU) 


The boundary conditions (55) require that 
f=) -3() for ¢=(0 and, |. Se tues at (GL) 
Now multiply eqn. (60) by F and integrate over the range of ¢ 
(ie. from 0 to 1). The left-hand side of the equation gives 
1 rl 
[ PD —ayF di= { FD*F dt—a? | Fd 
0 0 0 
1 
—[FDF}— [ {(DFP+aFhd. . . . . .- (62) 
~ 0 
By virtue of the boundary conditions on F' the integrated part vanishes 


and we are left with 


[ FD —ayF at=— [ (ore+ernat. fin. (63) 
0 0 
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Turning next to the right-hand side of eqn. (60), we have to consider 
2 Hh 
| ‘WF dt= i W(D2—a2)2W a—Q{ WD2W df... . (64) 
0 0 0 


Considering the first of the two integrals on the right-hand side of (64), 
we obtain after two successive integrations by parts 


| * W(D2—a2)[(D2—a2) Wat 
0 


ee | * WD[(D?—a2) W\dl—a? { * W(D?—a2)W dt 
0 0 

= {W D[(D?—a?) W|—DW[(D?—a?) W]} 

+] * D2W[(D?—a2) Wydl—a? | " W(D2—a®)W dt. . . (68) 
0 0 


The integrated terms vanish on account of the boundary condition W=0 
and either of (54). Thus, 


| W(De—a8)?W d= | [(D*—a2)WPdl. . . . (66) 
On the other hand 
| WDW dt=[WDW}— IF (DW)dt—— | (DW)'dt. . (67) 
Combining (64), (66) and (67), we have 
[we ac= | {(D*—ayWE+QDW)Yd. —. . . (68) 


The result of multiplying eqn. (60) by F and integrating over ¢ is, 
therefore, 


i (DF)--a2 Pal 


R= 


Te (say). ee (69) 


~ @2 


a | (Da) WP+-QDW)}ae 


This formula expresses R as the ratio of two positive definite integrals. 

Consider now the effect on R of an arbitrary variation SW in W 
compatible only with the boundary conditions on W. We have, to the 
first order 


I ds 
R= ar (8h 7h) = ap @I, (6J,—a*®RoI,), 2 eae 
where 5/, and 83/, denote the corresponding variations in J , and I,: 
51,=2 if {(DF)(D8F)--a? FSF} dt 
and » (iy 
‘1 
51,=2 | {((D*—a2) WI[(D*—a2)8W}+- Q(DW)(DBW) ak, 


Inhibition of Convection by a Magnetic Field 511 


The expressions for 5/, and SJ, can be reduced by integration by parts 
and making use of the boundary conditions 


. SW=8F=0 for £=0 and 1 


and 
D?sW=0 _ ona free surface, and ete) 
DsW =0 _ ona rigid surface. 
Thus 
1 1 
61, =2[6F DFY—2 | OF (D?—a?)F df= —2| 6F(D?—a?)F df. (73) 
0 0 
and 


51.=2{DW|(D2—a2)8W]—WD[(D2—a2)3W]--QW DSW}! 
+2 W {Da QD™pW}at 
—9 fe Wsredt. : 
Bil : Seer Ss ten ole Tivet, va FE) 
Now combining egns. (70), (73) and (74), we obtain 
29 
3R=— ar [ ard ar -+eRW}ae, Nes ee 75) 


Hence, to the first order 6k=0 for all small arbitrary variations in W 
satisfying the boundary conditions provided 


(D2 a2) Fone RW. == 01h otters, 21 (76) 


ie. if the differential equation governing W is satisfied. It is evident 
that the converse of this proposition is also true. Further, it follows 
from (75) that the true solution of the problem leads to a minimum value 
for R when evaluated according to formula (69). This last fact enables 
us to formulate the following variational procedure of solving eqn. (51) 
(for any assigned a?) and satisfying the boundary conditions of the 
problem. . 

Assume for F an expression involving one or more parameters, A,, 
and which vanishes at =0 and 1. With the chosen form of F determine 
W as a solution of the equation 


fe eae DW eee cee, (7%) 


and satisfying the boundary conditions on W at ¢=0 and 1. Then 
evaluate R according to formula (69) and minimize it with respect to the 
parameters A,. In this way, we shall obtain the ‘ best’ value of R for 
the chosen form of F. We shall see in the following section that even 
with the simplest trial functions we reach quite high accuracy in the 
deduced values of PR. 
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§6. THE DETERMINATION OF Rk, WHEN H AND g ACT IN THE SAME 
DIRECTION AND FOR THE CASES (i) WHEN BOTH BOUNDING 
SURFACES ARE FREE, (ii) WHEN BOTH BOUNDING SURFACES ARE 
Ricip AND (iii) WHEN ONE BounpING SURFACE IS FREE AND THE 
OTHER IS Rigip 

(i) Both bounding surfaces free —In this case the solution of eqn. (51) 

satisfying the boundary conditions (cf. eqn. (56)) 


W=D?W=D4W=0 for ¢=0 and 1, pn a 
can be explicitly written down. Thus 
W=constant sin nz, -‘eigSeeeaee (79) 


where n is an integer clearly satisfies the boundary conditions (78) and 
represents a solution of eqn. (51) provided 


R= (nx? a2) | (n?n? + a2)? + Qnen, eo GN 


_ That the solution obtained in this manner by inspection is in fact the 
unique solution of the problem can be confirmed by verifying that 
formula (69) also leads to the result (80), if the choice 


F=Asinn7l (A=constant), J Se Se 


be made for the trial function and the corresponding solution of eqn. (77) 
which satisfies the boundary conditions (75), namely 


A sin n7l 
w= (a2) + Onn > . . . ° . . (82) 
be used in the evaluation of R. 
Letting a?=n?7x, we can rewrite eqn. (80) in the form 


R=ntnt —* [ota at 


a= (83) 


From this equation it follows that for a given 2 (i.e. a?) instability first 
sets in for the lowest mode n=1. The corresponding expression for R is 


lt+2 ¢ Q 
R=nt —— [1 +2)?+Q,1] (2.5). - « » (84) 
As a function of x, R given by the foregoing formula attains its minimum 
when 


2a3 + 32 =O) ae leo i heen 


With x determined as a solution of this cubic equation, (84) will give the 
critical Rayleigh number, R,, at which convection through instability 
will first set in. The values of R, determined in this fashion for various 
values of @ are given in table 1. The inhibiting effect of a magnetic 
field on the onset of convection is apparent from this table. 

For @, sufficiently large the root of eqn. (85) tends to 


Cin >(Q/2)*8 ; oo Ae me se ek) 
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and the corresponding limiting expression for R, is 


Resta] ep tth e4 e s -(87) 


It is found that this last formula for R, gives sufficient accuracy for 
Q>1000. 


Table 1. Critical Rayleigh Numbers for the Case when both Bounding 
Surfaces are Free 


Q a R, R.|Ro Q a R. R.|R, 


0 | 2-233 657-511 1-0000 6500 7-952 80343-6 | 122-19 
5 | 2-432 796-573 1-2115 7000 8-059 86034:0 | 130-85 
10 | 2-590 923-070 1-4039 7500 8-159 91705:7 | 139-47 
20 | 2-826 1154-19 1-7554 8000 8-253 97360-1 | 148-07 
50 | 3-270 1762-04 2-6799 8500 8-343 | 102998-0 | 156-65 
100 | 3-702 | 2653-71 4-0360 9000 8-429 | 108623-0 | 165-20 
150 | 3-990 | 3475-67 5:2861 9500 | 8-510 | 114234-0 | 173-74 
200 | 4-210 | 4258-49 6-4767 | 10000 8-588 | 119832-0 | 182-25 
300 | 4:543 | 5752-65 8-7491 | 10500 8-663 | 125419-0 | 190-75 
400 | 4-794 | 7185-94 10-929 11000 8-735 | 130995-0 | 199-23 
500 | °4:998 | 8578-88 13-048 11500 8-804 | 136560-0 | 207-69 
600 | 5-171 9942-40 15-121 12000 8-870 | 142116-0 | 216-14 
700 | 5-321 | 11283-2 17-160 13000 | 8-997 | 153202-0 | 233-00 
800 | 5-455 | 12605-6 19-172 14000 9-116 | 164255-0 | 249-81 
1000 | 5-684 | 15207-0 23-128 15000 | 9-227 | 175279-0 | 266-58 
} 1500 | 6-123 | 21535-2 32-753 16000 9-333 | 186276-0 | 283-30 
| 2000 | 6-453 | 27699-9 42-128 17000 9-433 | 197249-0 | 299-99 
2500 | 6-720 | 33756-5 51-340 18000 9-528 | 208199-0 | 316-65 
3000 | 6-945 | 39734-2 60-431 19000 | 9-619 | 219129-0 | 333-27 
| 3500 | 7-140 |45650-6 69-429 20000 | 9-706 | 230038-0 | 349-86 
| 4000 | 7-313 |51517-8 78-353 25000 | 10-09 284341-0 | 432-45 
4500 | 7:442 |57344-6 87-215 30000 | 10-42 338308-0 | 514-53 
5000 | 7-585 | 63135-9 96-022 35000 | 10-70 392013-0 | 596-21 
| 5500 | 7-717 | 68897-3 104-78 40000 | 10-95 445507:0 | 677-56 
6000 | 7-839 | 74632-1 113-51 


(ii) Both bounding surfaces rigid—wWe shall obtain the solution for 
this case by an application of the variational procedure described in § 5. 
But first it is convenient to translate the ¢-axis so that the bounding 
surfaces are at £=+4. With the bounding surfaces at C[—-+-3 it follows 
from the identity of the boundary conditions to be satisfied at €=-+4 
that the solutions for F and W should be either even or odd functions 
of £; and it can be readily verified that the even solutions give lower 
values for R. We must therefore choose the even solutions. 

As a simple trial function for F which satisfies the boundary condition, 
F=0 at €=+4, we shall assume 


F=A cos wf+1+ cos2al (—4<0<+4), - « . - (88) 
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where 4 is a variational parameter. With F given by (88), the equation 
governing W is 


[(D?—a®)?— QD?] W=1-+-cos 27+ A cos mC; i ae Bas 
and the solution of this equation appropriate to the problem ‘on hand is 
W=yo+y-2 cos 2n6-+ Ay, cos f+ B cosh g,6+C cosh got, . (90) 
where B and C are constants of integration, 
H=HV(Q+4e*)+ VQl, g=sV(O+407)—VQ], - - (91) 
and 
(92) 


The boundary conditions to be satisfied by W are (cf. eqns. (53) and (54)) 
We DW =O Moret =a. Tee Lhe ee eee 
These conditions require that 


B cosh 3q,+C cosh $¢2=—(vo—y2) 


ray 


and ae eee 
Bq, sinh 4¢,+Cq,_ sinh $q¢,.=Az7y). 


The solutions of these equations can be expressed in the forms : 


B= + 0[Ary,+4o(yo—y2) tanh 349] sech 3q;, (95) 
5 
C=—8[Ary,+9i(yo—yz2) tanh 39,] sech 3q2, 
where 
; 5=(q, tanh 4¢,—q, tanh }q.)"7. . . . . (96) 


With F and W given by eqns. (88), (90) and (95) we find that 


+1/2 39 472-1. 32 
bse [(DF)?+a2F?]df=4(7?+-a?) | ae+ At a - (97) 


and 


+1/2 47728 ( (q 2_g 
WF dt=} E 142") ae | az 
fee Sin (7? -++-q,")(7?-+q0" 


‘ tanh 1 t nh 1 
sae E (38yo +471 +Y2)—278(yo—y2) { 71 COD ad ste ee 


+o? ge 
tanh $q tanh 4q 
as 24 i es 291 
rae { J2(47° +42") eee || pa | (vot Aya) 
(41° —2")(47?-+-9"-+- 45") 
a(47?+-9,°)(402 +452) (98) 


The resulting expression for R is, therefore, of the form 


4-1/2 
DF)2--a2F2 
Be J _lDP er Wo f+gA+ha® 
~ l+mA+nA?2° 


— 87?(yo—yo)d tanh $q, tanh 292 


+2 (99) 
ae WF dé 


—1/2 
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As a function of A, R attains its minimum when A is a root of the equation 
(gl—mf )-+2(hl—fn) A + (hm—gn)A?=0. ~'. . (100) 
With A determined by this equation, eqn. (99) determines the ‘ best ’ 
value of R for any given a. By repeating the calculations for various 
assigned values of a we can determine the value of R,. However, in 
practise it was found more convenient to carry out these calculations in 
two stages in the manner described below. 
First, we may observe that if we had assumed for F the trial function 


F=A cos xé, See ees se) ewe! (101) 
we should have found for R the formula 


_ a 4am? +.02)?-+-299)] 

a1 — 41770 (44?— 42?) /(7? +-47)(7? +95") | ; 
which follows from eqns. (97) to (99) by letting A—oo. It appears that 
this formula which has no variational parameter in it already determines 
& to an accuracy better than 1%. For in the field-free case when Q—0 
and ¢,;=q.—a the limiting form of (102), namely 

(7? a?)8 
~ @{1—167?a cosh? da/(a-+ sinh a)(7?+-a?)?]’ 
gives &,=1715-1 for a=3-13; and this value should be compared with 
&,=1707-8 determined by Pellew and Southwell (1940) by an exact 
solution of the problem. The procedure then was the following : 

For an assigned @ the value of a which gave the minimum value for R 
evaluated according to eqn. (102) was first determined. With this value 
of a a ‘second approximation’ was then carried out with the aid of 
eqns. (97) to (99) by making the variation with respect to A. The 
results of such calculations are summarized in table 2. The extent of 
the agreement between the results of the first and the second approxima- 
tions may be considered as a measure of the accuracy of the determination 


Table 2. Critical Rayleigh Numbers and Related Constants for the Case 
when both Bounding Surfaces are Rigid 


(102) 


R (103) 


First Second 
Q@ g “a Approximation | Approximation Hal Feo 
0 3-13 1715-1 1707°8 1-0000 
10 3°25 10-506 1953-7 1945-9 1-1394 
50 3°68 10-963 2811-4 2802-1 1-6408 
100 4-00 11-687 3767-6 3757-4 2-2001 
200 4-45 12-910 5499-9 5488-6 3-2138 
500 5-16 FSO Pe 101220 10110-0 5:-9199 
1000 5-80 19-261 17116-0 17103-0 10-015 
2000 6-55 | 23-830 30139-0 30125-0 17-640 
4000 7-40 29-896 54711-0 54697-0 32-028 
6000 7-94 34-500 78405:0 78391-0 45-902 
8000 8-34 38-052 101622-0 101606-0 59-495 
10000 8-66 40-462 124523-0 124509-0 72-906 


pO 
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of R, by the variational method. In this connection it should be noted 
that an uncertainty in a of 0-1 in the neighbourhood where #, attains 
its minimum does not affect the deduced values of R by more than 
one part in 10°; and often, it is much less. 

Finally, we may remark that the asymptotic behaviour of F for large Q 
is the same as in case (i). 

(iii) One bounding surface free and the other rigid.—In this case the 
conditions to be satisfied on the two bounding surfaces are different. 
But as has already been remarked by Low (1929) and Pellew and 
Southwell (1940) the solution for this case can be reduced to that of 
case (ii) by considering the odd (instead of the even) solutions. For, it 
is evident that solutions for F and W which are odd and satisfy the 
boundary conditions of case (ii) at ¢=+4, vanish at ¢—0 and satisfy 
here the boundary conditions appropriate to a free surface. Consequently, 
an odd solution for case (ii) applicable to a cell depth d also provides a 
solution for case (iii) and applicable to a cell depth d’=4d and a Rayleigh 
number which is 16 times smaller. 

In accordance with the foregoing remarks we shall assume for F the 
trial function 
F=A sin 27f+sin 7f+ sin 3r€. hoe eye . (104) 
The corresponding solution for W is (ef. eqn. (90)) 
W=y, sin 7f+y, sin 37€-+-Ay, sin 27f+-B sinh g,€+C sinh q,., 
ey (105) 


where y,, q; and q, have the same meanings as in eqns. (91) and (92). 
And the integration constants B and C determined by the boundary 
conditions (93) are: 


B= +8'[2ry.4 +9o(y1—ys) coth }qy] cosech 291) 

C= —8'[2ry.4 +9;(y1—ys) coth $q,] cosech 445, 
where 6’=[q, coth 4q¢,—q. coth4g,]-?. . . ... (107) 
With F and W given by the foregoing equations we find that 


me 2(57?-- a?) 


(106) 


+1/2 
i [DFP +P db=1(4n?-+0°) 


and 


e+1/2 pos 2 2 
| WF dt=1y, E seve Be 
—1/2 (47°-+-q,")(4?+-95”) 


4 q, coth $ 4 
Be te Yee q, coth 39, — Jz coth 3g 
| On t 8y2+ 873) —478"(y, rs) “4g ercornr 
th 4g, th 4, 
39738) 2 ; Yo COUN 2do rs Jy COUN $94 
LMP) GH Onetgays lA + | Matys) 


2 ps 2 2 2 
—- 1672(y py; )d’ coth 4q eoth 4 $ "Oem a0 tae tae) 
ae OS BELG VE GEO GO Pay) | 


(109) 
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As we have already stated R determined by formula (99) and eqns. (108) 
and (109) is 16 times the value of R for the case under consideration. 

As in case (ii) the determination of R, in accordance with eqns. (99), 
(108) and (109) was carried out in two stages : an approximation with no 
variational parameter and corresponding to a trial function 


F=A sin 2rC, hs GMS ee Se ei) 


was first carried out. The formula for R in this first approximation 
(which can be obtained by letting A—oo in eqns. (108) and (109)) is 
fa (4n2-+-a2)[ (47? +22)? +472] 
— 16a*[1— 1618" (q?—957)/(40?-+9,7)(47? + 90")]° 


(111) 


Table 3. Critical Rayleigh Numbers and Related Constants for the Case 
when One of the Bounding Surfaces is Free and the Other is Rigid 


First Second . 
@ 2 as Approximation | Approximation Fie| Ro 
0 5:36 4-89 1112-7 1100-75 1-0000 
10 5:50 4-96 1179-4 1167-2 1-0604 
50 5:94 5:25 1428-3 1415-5 1-2859 
100 6-34 5:59 1712-7 1699-4 1-5439 
200 6-9 6:19 2231-3 2217-6 2-0146 
500 8-0 7-54 3600-2 3586-1 3:2579 
1000 9-0 9-21 5627°5 5613-3 50995 
2000 10-2 11-6 9318-7 9304-5 8-4529 
4000 11-5 14-9 16133-0 16119-0 14-644 
6000 12-4 17-3 22606-0 22592-0 20-524 
8000 13-0 19-4 28893-0 28879-0 26-236 
10000 13-5 21-2 35058-0 35044-0 31-836 
20000 15-3 27-9 64861-0 64847-0 88-912 
40000 17-3 37:0 122155-0 122140-0 110-96 


The value of a for which this formula gave a minimum value for & was 
first determined. Then with this value of a a second approximation was 
carried out with the aid of eqns. (108) and (109) by making the variation 
with respect to A. The results of such calculations are summarized in 
table 3. 

It may be noted here that in the field-free case when Q@=0 and 
91 =%.=« the limiting form of (111), namely, 
(4a? +a?) 


Le 16a?[1—647°a sinh? 4a/(sinh a—a)(47?+-a?)?] 


(112) 
gives R,=1112-7 for a=5-36; and this value should be compared with 
R,=1100-65 determined by Pellew and Southwell (1940) by an exact 
solution of the problem. On the other hand a second approximation 
carried out with the corresponding limiting forms of eqns. (108) and (109) 
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gave the value 1100-75. As in case (ii) we may, therefore, expect that 
the tabulated values of R, in the second approximation are correct to 
at least one part in 10°. 

Finally, it may be remarked that for sufficiently large Q the values of 
R, for this case approach a quarter of that for case (ii) (or (i)). 

The results of the calculations for all three cases are illustrated in fig. 1. 
It is apparent from this figure, as well as from tables 1 to 3, that the 
existence of a magnetic field will have, depending on the value of @ and 
the boundary conditions, an inhibiting effect on the onset of convection ; 
and this inhibiting effect will be the greater the greater H and (or) o 
are (is). In a general way the reason for this inhibiting effect is clear. 


Fig. 1 : 


log Re 


4:0 


20 


a) 2:0 3:0.) 4:0 
log Q 
The variation of the critical Rayleigh number R,, for the onset of instability as a 
function of @ for the three cases (i) both bounding surfaces free cae 
marked 1), (ii) both bounding surfaces rigid (curve marked IT) and (iii) on 
bounding surface free and the other rigid (curve marked III). The Soe 
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When the field is strong (or the conductivity high) motions at right 
angles to H will become increasingly ~ difficult > and this will prevent 
an ‘easy’ closing in of the stream lines required for convection. This 
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also explains why when instability finally sets in, it does so for a value of a 
which is increasingly large: for, a large a means that the cells are 
elongated in the direction of H (cf. eqns. (46) and (49)) and motions in 
this direction are not hindered. 

In order to see how effective the magnetic field will be in practical 
cases in inhibiting the onset of convection we shall consider mercury at 
room temperature. We have 


G11 xX 100) ppeal “7 x 10-2 
and 


p= Hod? 


where it may be recalled that H is measured in gauss. For 


=6% 10-4 Heda, ee ns. (118) 


d=1, H=103, Q=600: ee te eee (114) 
and the critical Rayleigh numbers for the three cases considered are 
Ki= 9 942, &,/Ro=15-12, in case (i), 
= 11 6605 4 / ho 6-80, in case (11), ee ell) 
;—= £018, R,/hj= 3-65, im case (i). 
It would appear from these values that the predicted effect should be 
easily detectable in the laboratory. In §9 we shall indicate certain 


astrophysical and geophysical connections in which this effect may 
play a role. 


§7. THE CASE WHEN H AND g AcT IN DIFFERENT DIRECTIONS 
We shall now briefly consider the case when H and g act in different 
directions. 
Let H act in a direction inclined at an angle # to the vertical, i.e. to 
the direction in which g acts. Also, let the direction of the x-axis be so 
chosen that H lies in the xz-plane. Then 


A= (On0, by andet=.(sind. 0. cos), . 1. . .. (116) 
and eqns. (38) and (39) become 
Vip y é& 3) pHa 
Vv 


Vv 


oa? | ey? 


fo] = 0 2 lod 
(cos 5. + sin 95°)" w y Lt?) 


and 


ioe a 
7°y— +—— | —~ ——— 


0 ; a \2 
KV Ox? oy? + (cos 95: sin 9 ) V0, ~(118) 


ox 
If we seek solutions of eqns. (117) and (118) which are independent. 
of x and are of the forms 


w= W(2) cos + and 6=O(2) cos, vn eee oe 119) 


2N2 
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then eqns. (117) and (118) reduce to eqns. (50) and (51) with the only 
change that in the definition of Q (eqns. (52)) H cos d replaces HH, 
Consequently, if we restrict ourselves to an onset of convection in the 
form of infinitely extended ‘rolls’ (instead of cells; in this connection 
see Chandra 1938) the discussion of the criterion for stability can be 
carried out exactly as in §§4 to 6; and the critical Rayleigh numbers 
can be deduced from those given in tables 1 to 3 by interpreting H to 
mean the component of H in the direction of g. The question whether 
lower Rayleigh numbers can be reached by considering patterns of 
motion more complicated than rolls will be left open for the present 
though it would appear that on the whole the possibilities for this 
happening are rather remote. 


§ 8. AN EXAMINATION OF THE PRINCIPLE OF THE EXCHANGE OF 
STABILITIES 


As we have already pointed out, the description of the situation in 
marginal stability in terms of the relevant equations in which o/dt has 
been set equal to zero is valid only if the principle of the exchange of 
stabilities can be justified. Now this principle supposes that in the 
neighbourhood of the state in marginal stability the effects of a free 
disturbance of the steady state will vary like e® where w is a real number 
so that o/ot—w=0 does define the marginal state. The principle will 
also be valid if it can be shown that even if w is complex, the limit of the 
situations with a positive real part for w is one in which both the real 
and the imaginary parts tend to zero simultaneously. But the principle 
will not be valid if when the real part of w tends to zero from the positive 
side (i.e. &(w)—>-+0) the imaginary part of w tends to a finite limit. In this 
latter case when instability develops, the ensuing motions will be of 
the nature of oscillations with increasing amplitude, i.e. we shall have 
‘ overstability ’ in the sense of Eddington. 

For the classical Rayleigh—Jeffreys problem, Pellew and Southwell 
(1940) have proved by a very general argument that for all adverse 
temperature gradients (i.e. B=—| |) the time dependent equations admit 
solutions varying like e” with w necessarily a real number. The 
application of the principle of the exchange of stabilities is therefore a 
valid procedure for that problem. On the other hand it will appear 
that the method of Pellew and Southwell is not quite strong enough to 
establish the principle for the problem on hand. We shall therefore 
examine a particular case in some detail and then see how far the general 
argument enables us to extend the deductions drawn from the particular 
case. 

We start then from the general time dependent eqns. (28) to (30). 
Kliminating w from the first of these equations by operating with V? 
and making use of eqn. (26), we obtain 

0 te: pH 020. 


0 
= \ 2) ee ee roe ee 
or W220 a= dia, an, Vh+y (v9 XA; sat) EVV S10 ie mea Doe 
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If all quantities vary like e” the foregoing equation and eqns. (29) 
and (30) reduce to 


H 0 076 
ee NO\ a et 
(w—vV?)V2w= ECA Ae Why (V8 Vii aap AB ), 


(w—xV?)0=—Bw, ma Lo) 
aie 

t Ox, 

In the further ieee of this problem we shall restrict ourselves to 


the case when both H and g act in the direction of the vertical (OMO REI): 
Then eqns. (121) become 


(w—nV?)h=H 7, 


ne; Gath. ofa) 
<A ee 
Anion” vinty (53+ 5p): 
ae eee (to) 
(w—nV2)h= ee 


By operating the first of these equations by (w—«V?)(w—7V?2) we can 
eliminate h and 6 and obtain an equation only for w. Thus, 


2 92 
(w—KV?)(w—7V?)(w—rV?2) V2wW= eb 


Apo Oe (w—KV?2)V2w 


oon 0? 
We shall reduce this equation still further by assuming, as in § 4, that 
the variables x and y can be separated from z and that w is of the form 
w=f(a, y)W(z) meray Meo ets. ca (124) 


where f satisfies eqn. (44). For w of this form, eqn. (123) can be 
reduced to 


[Q—K(D?*—a*)][Q—n(D*—a2)][Q—v(D*—a)|(D?— a2) W 
oe "| Q—K( D2—a2)]D(.D?—a2) W—a? a B | y[Q—n(D®—a2)] W, . (125) 


ee Oeawdrat = (inane —=Ol0l. “0 «.--..27. «= (126) 


A solution of eqn. (125) must be sought which satisfies the various 
boundary conditions of the problem. By letting 


We_constantistn gre we. : 4. «+ x (127) 


we can satisfy the boundary conditions appropriate to both bounding 
‘surfaces being free (i.e. 0= W=D?W=0 for ¢=0 and 1) and further have 
the normal component of the gradient of h vanish on the bounding 
surfaces. A straightforward discussion of eqn. (125) for any other set 
of boundary conditions is likely to be troublesome. It will not be 
attempted in this paper. 


t 
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For W given by (127) eqn. (125) becomes 
[Q+(2?-+-a?)] [Q+-y(7?+4?)] [Q+(7?--a*)] 


pd? : a2 : ; 128 
a 7 O+«(72+a?)]—|B ly dt [Q+9(7?+0?)]=0. (128) 
dit pg a?+o 
Letting 
Q d? a 
i haps Fae and y= 4 tl, ar ee 4 5 a A (129) 


we can rewrite eqn. (128) in the form 


(qtKy)(q+ny)(q+ry) +O inv(q+Ky) 


ii ee 
—Ryw(q) a), ites is oe (130) 
where (cf. eqns. (52)) 
Q H?o d* Rk |Bly@ 
=3=—s = 5[="— jn ee 
a= =p? aes and k, ms pert (131) 


Equation (130) is a cubic equation for g. Writing it in the form 


g++ Bq? Cg D=0; 5. 5s eee 
we have 
D=nnlyY+Qy—Ryy—V}, | 
Viz : 
(eae) nye Ra a eo cath ane 
B=ys1, §;=K+v+7 and sg=Kv-+vn-+y7k. | 


From eqn. (132) it follows that if D<0, it admits a positive real root. 
Consequently, if D<0 and we let it tend to zero from the negative side 
we shall have g—0 through positive real values. In other words, in this 
case we shall approach marginal stability as the limit of solutions which 
diverge as e” with w real and positive; D<0 is therefore a sufficient 
condition for the principle of the exchange of stabilities to be valid. 
Therefore, let 


D>0. 50 kT eee ee rie, ee ae oem 


Equation (132) then allows a negative real root and the corresponding 
solutions lead to an exponential damping of an initial disturbance ; 
these solutions (as well as their limit D->+0) are of no interest to us. 
We have now to examine whether under these conditions the complex 
roots of (132) (if such exist) can have a positive real part. (For a similar 
discussion of a cubic equation see Jeans (1929, p. 120).) If p/2 denotes the 
real part of a complex root of eqn. (132), then it satisfies the equation 


p®-+2Bp?+ (B2+C)p+(BC—D)=0. . . . . (135) 
If now BC—D<0, cient Sahel, acts Wee ks ay Oe 
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eqn. (135) will admit a positive real root and by letting BO—D-—-0 
through negative values, we shall have the case when a complex root of 
eqn. (132), having a positive real part, tends to a purely imaginary limit. 
Thus, if 


D>0 and BC—D<0, ee tae e137) 


the principle of the exchange of stabilities cannot be applied and we will 
have a case of over-stability. 
According to eqns. (133), the conditions (137) are equivalent to 


Ry(y—1)<y?+Qiy ers Lo) 
and 
co {yP+Q,y—R,(y—1)}>y3s,8.+Qinvs,y—RKvs(y—1). . (189) 


Rearranging this last inequality, we have 


1 
Se ah) heap) {y?(818.—Kvn) + Qi nr(ntv)y}. . . (140) 


The inequalities (138) and (140) can now be combined in the form 


P+Qy>Riy—1) > (yee) + Qin typ}. (141) 


KV 


Since y is necessarily positive (cf. eqn. (129)), we can rearrange the 
outer inequalities of (141) to give 


Q, {1— vo Se ee Mrtuteletaltze]. . (142) 
This last inequality is a necessary condition for instability to arise through 
oscillations of increasing amplitude; it is therefore a condition for the 
principle of the exchange of stabilities not to be valid. 

According to (142) a sufficient condition for the principle of the exchange 
of stabilities to be valid is 


re) or 7>k. peer et tae wica (14a) 


(In § 10 we shall show that this condition is also necessary.) It would 
appear that for fluids under terrestial conditions this last inequality, 
7 >«, will be satisfied with a large factor of safety. Thus, for mercury 
at room temperature 


n=7-5 X 108 cm?/sec and «=4:710-? cm?/sec. . (144) 


Accordingly, in these cases we may safely apply the principle of the 
exchange of stabilities. However, under astrophysical conditions 
(see § 9) <x and we must go back to the general inequality (142); we 
shall return to this in § 9. 
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If condition (142) should be met and we have a case of over-stability, 
then the frequency of oscillation w’ at marginal stability will be given by 


7 KVY) 1/2 
oF [y8+Qiy—R,( y—1))| 


YS, 


7 Yai) 12 a 
= yt OR | : ss tag aw Len 


and the critical Rayleigh number for the onset of convection can be 

obtained by minimizing (with respect to y) the expression 

gcic2— KVR P+Oy 
Kv(kK+v) y—l 

Q n(q+r) 


7 8)89—Kvn 


R= Sage nee CAG 


where aN (147) 

The foregoing discussion applies to solutions of eqns. (122) which 
satisfy only a particular set of boundary conditions. For solutions 
satisfying more general boundary conditions we might treat eqns. (122) 
in the manner of Pellew and Southwell (1940, §§ 16 to 18). Thus if w is 
complex, let w*, 6* and h* belong to the complex conjugate, w*, of w. 
Then (cf. eqns. (122)) 


we ye BES rapacaie | Clea 
(w* —pV2)V20 a eae? (Sp ip 


Multiplying this equation crosswise by the second of the eqns. (122) 
we obtain 


(148) 


* 29% 
Bw(wo* —vV2) V2w* — ; a we V2h*+-y(w—KV2)0 (Ge +r) =0. (149) 
Now we integrate this Lae over a unit cell. The various volume 
integrals can be reduced by transforming them in accordance with 
Green’s theorem and making use of the proper boundary conditions on 
the cell walls and the bounding surfaces. The integrals arising from 
the first and the last terms on the left-hand side of (149) have already 
been reduced by Pellew and Southwell. It remains to consider 


I] we Ve dedydz. . . . . . .. (160) 


Integrating this by parts (with respect to z) and remembering that w 
must in all cases vanish on the bounding surfaces, we have 


—[ [| 2 ven de dy ae see sae CT 


Replacing dw/dz by (w—nV2)h/H (ef. the last of the eqns. (122)) and 
transforming the integral by Green’s theorem, we obtain 


1 ip é Sor 
All| {o| grad h |?-+-n| V2h |?} da dy dz— al [a grad h* dS. (152) 
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The surface integral in (152) vanishes since on the cell walls (being 
surfaces of symmetry) oh*/On vanishes and on the bounding surfaces 
(¢=0 and 1) either / or dh/dz must vanish. We are therefore left with 


(ero. 1 
[| fos ven da: dy dz= lwlgtyl,), . . . . (158) 
where J; and /, denote the positive definite integrals 
i {} | grad h |? dx dy dz and I,—= {| [ [V2 pe de dy de. (154) 


The result of integrating eqn. (149) over a unit cell is therefore (cf. Pellew 
and Southwell 1940, eqn. (29)) 
ee up 
Berl + F Tetyolat Polat Kyl + oe Ig=0, . . (185) 
where 
ed | [grad wh dx dy dz, 1, || [ive [2 dx dy dz, 
a ihe aa 


io = {|{ [0 |? da dydz and I,= = ||| lerad op da dy dz 


are all positive definite integrals. 
Separating the real and the imaginary parts of eqn. (155), we have 


i iI 
R(o) {71,48 (1+ £2) +6 (tHE!) tor—0 . (197 


0 


and Ho) | —B (1, E*) +71} =0. eee re LOS) 
Po 


From eqn. (158) it follows that %(w) is necessarily negative if P>0; 
in other words when the temperature gradient is not adverse all small 
disturbances will decay. For (w) to be positive it is a necessary 
condition that 8 be negative. Let B=—|f|. Then eqns. (157) and (158) 
give 


| B |(vLo+wnle/47po)—KyLy 
SO) a Mt ee Oo 
oy yl;—| B (1 +4-15/477p0) Oy 
pls ) | 
d I,-— 1 ely ante en LOO 
an Ho)| 121( 1 4irp, +yty (160) 
From egn. (160) we conclude that 
a te 2 
evther J(w)=O0 or 18 1(,— fs) eye On meant (LOL) 


In the former case we recover the principle of the exchange of stabilities 
and we may readily verify that eqn. (159) then reduces to the formula (69) 
for evaluating R in the case of marginal stability. In the latter case, we 
may have over-stability; and a necessary condition for this to be 
possible is 


eg 
—— ge ee eee (162 
aes Tees, (162) 
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or, explicitly (cf. eqns. (154) and (156)) 
ie \|] | grad h |? da dy dz> | || [grad w[?dadydz. . (163) 
ye it 
It would appear that this last condition can in general be satisfied only 
if the magnetic energy associated with the disturbance exceeds the 
kinetic energy. Under most conditions we should not expect this to 
happen: for, even in a magneto-hydrodynamic wave the two energies 
are only equal. We conclude that it is highly probable that the principle 
of the exchange of stabilities will be valid under most terrestial 
conditions; but it is by no means certain that this will be the case under 
astrophysical conditions. A further examination of the inequality (142) 
in § 9 in the light of (163) confirms this conclusion. 


$9. Tur INHIBITION OF CONVECTION BY A MAGNETIC FIELD 
UNDER ASTROPHYSICAL CONDITIONS 


On the grounds that motions ensuing on the onset of convection must 
be characterized by closed stream lines and that such motions would be 
hindered by a magnetic field, it has been suggested (cf. Biermann 1941 
and Walén 1949*) that the presence of a magnetic field might suppress 
convection altogether. Indeed, one has sought to explain the appearance 
of sunspots in terms of such an inhibiting effect of the magnetic fields 
associated with them. ‘The theory developed in the preceding sections 
enables us to examine to what extent these suggestions are sound. 

Consider, then, the conditions in the so-called hydrogen-convection 
zone in the solar atmosphere where 


pai2X10-' g/em*, T=6 0007 kee ee 


and the most abundant constituent of the atmosphere, namely, hydrogen 
is getting about half ionized. Osterbrock} (1952) estimates that under 
these conditions 


x=10'8 cm?/sec, v=5 x 108 em?/sec and 77 X 108 cm?/sec. . (165) 


It will be noticed that contrary to terrestial conditions (cf. eqns. (113) 
and (144)) the largest of these coefficients, by far, is the thermometric 
conductivity, «. The reason for this is that under astrophysical conditions 
the transport of heat by radiation is an extremely efficient process. 
Indeed, it would appear that 


K>n>v, viel eee cle) a el eeL GG) 


is characteristic of astrophysical conditions. We shall now examine in 
some detail the extent and manner of the inhibitory influence of magnetic 
fields on convection when (166) holds. 


ng op eee 
* Tam indebted to Professor T. G. Cowling for drawing my attention to these 

references. aoe 

_ + 1am grateful to Mr. D. Osterbrock for allowing me to use the results of his 

investigation before publication. 


Inhibition of Convection by a Magnetic Field 527 


The first question to settle is of course the applicability of the principle 
of the exchange of stabilities. For this purpose we go back to the 
fundamental inequality (142). It will be recalled that this inequality 
was derived for boundary conditions appropriate to both bounding 
surfaces being free; it would not seem that the use of these boundary 
conditions for the astrophysical problem is artificial. 

Now when (166) holds, the inequality (142) can be simplified to 


Cy nee 2) 2.2 (167) 


Suppose that this condition is satisfied and that over-stability arises. 
As we have already pointed out, under these circumstances the critical 
Rayleigh number at which instability will first develop can be found by 
minimizing the expression on the right-hand side of eqn. (146). And since 
($8.—Kvy) ~k2y, when (166) holds, we can simplify eqns. (146) and (147) 
to the forms (cf. eqns. (131)) 


_ ny +Vy Pele 
tke eee ee Be where W'=41 3: eames ere (1O8) 


For Q, a multiple of y/v say An/v, Q’=An?/K?; we can, therefore, neglect 
Q’ in the expression for R, so long as A does not become comparable to 
k/n?, Consequently, for a considerable range of @, (satisfying (167)) 
we can write 


ar Tae a 169 
ie een (a.<*). eee eae ate A (LOO) 
The right-hand side of this expression has a minimum for y=3/2. 
Therefore, under the conditions specified, the critical Rayleigh number 
for the onset of over-stability is 


27 rie 


R.=7*Ry, min, =? (170) 


2 
Vv 


and this is independent of the strength of the prevailing magnetic field. 
Returning to the more general inequality (141), we observe that under 
marginal conditions this must reduce to 


PLOY > RY) ers) + Om ory}- AT) 


_We have seen that when (166) holds and Q,<x?/nv marginal stability 
obtains for y=3/2 and the term in Q, on the right-hand side of the 
equality in (171) can be neglected. Consequently, for 


1 <n/v<Q,<x*/nv, Ros hy, es Alle i eT a.) 
(171) becomes 
279 : 
$Q,>¢4;,, mines oe.) A ° ° . . . . é (173) 


and this is consistent with (170). 
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What the discussion in the preceding paragraph establishes is this: 
when (166) holds, instability in the form of oscillations of increasing 


amplitude can arise for 
9 
QS). OE 


and when it does develop by over-stability, the critical Rayleigh number 
remains practically constant and has the value 
27 
bp Tl aes. Me es Ce S) 
4 ov Vv 
so long as Oy <rtguil Ses ee Poe een a 


For Q,>+2/nv, R, increases in a manner which is identical with case (i) 
treated in §6; the results of table 1 can in fact be used quite generally 
if Q is interpreted as meaning 7? and FR, as the required Rayleigh number 
divided by (s,8.—kvn)/Kv(K-+y). 

On the other hand if Q,<9y/4v, the principle of the exchange of 
stabilities will be applicable and instability at the critical Rayleigh 
number will arise through ordinary cellular convection in accordance 
with the results of §6. Again, restricting ourselves to the case when 
both the bounding surfaces are free and the results for case (i) § 6 apply, 
we can, in the relevant range of Q,(~7/v>1), write (cf. eqn. (87)) 


R,=7'Q,. ya ei te), bee 


The value of Rk, given by (177) is less than that given by (175) so long as 


277 
0, <Q,*=—-. ge tay a ogee On ee 


For @,>@,* the reverse will be the case. We may therefore conclude 
that for Q,<Q,* instability will arise at the appropriate Rayleigh number 
by cellular convection while for Q,>Q,* it will arise through 
over-stability. 

We have seen that when the conditions (166) hold, the critical 
Rayleigh number for instability remains constant for a considerable 
range of Q, greater than @,*. But @, (ie. the magnetic field) influences 
the frequency of the oscillations with which over-stability sets in. Thus, 
according to eqn. (145), for Q@,>Q,* 


wo! = (Qn). vast 5) geen aE 


: : 
The meaning of this frequency becomes clearer when we express @, in 
terms of the velocity of the magneto-hydrodynamic wave, namely 


wH?\ 12 
i= (=): e . . . : ° ‘ (180) 


Thus (cf. eqns. (20) and (131)) 
pH dev a 
0, = =3. See poe ae fe tL 


— Arponvr nv 7 
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Substituting this last expression for Q, in (179) we find, 


7 
Q = d Vee ; 
2n/w’ is therefore the time required for a magneto-hydrodynamic wave 
to travel a distance equal to $d. Similarly, some insight into the meaning 
of Q,* (at which instability by convection changes into instability by 
oscillations of increasing amplitude) can be obtained by expressing it in 
terms of the velocity of the magneto-hydrodynamic wave. Thus, letting 


ee ae hee 


(182) 


* si 
7 ar epe i ee (183) 
we have 
4/27 7 
LoS ite at : 6 cee Go (184) 


We may now summarize the results of the foregoing discussion as 
follows: so long as the velocity of the magneto-hydrodynamic wave is 
less than V,,* given by (184), instability will’arise through cellular 
convection when the appropriate Rayleigh number derived on the 
principle of the exchange of stabilities is reached. But if V,,,*<V,,,<ak/d 
we will have over-stability when the Rayleigh number reaches the value 
274y/4v; and the frequency of oscillation at marginal stability will be 
essentially determined by the time required for the magneto-hydrodynamic 
wave to travel the thickness of the atmosphere. All these conclusions, 
it should be repeated, are valid only if the (166) holds. 

Returning to the numerical values (165) for the basic coefficients, we 
observe that 

eed Oc atone Ss <@L02 oe O21 084 ie a 2) (185) 


Therefore, under astrophysical conditions, the inhibition of the onset of 
instability by a magnetic field is not very pronounced. On the other hand 
since in all practical cases we may expect the instability to develop in 
the form of oscillations of increasing amplitude, the question whether in 
these cases we will have convection in the true sense is not certain. The 
possibility of explaining sunspots in the manner indicated at the beginning 
of this section must therefore be left open for the present. 


§10. THE MANNER OF THE ONSET OF INSTABILITY IN THE 
GENERAL CASE 


In §§ 8 and 9 we have seen how in virtue of the inequalities 7 >« and . 
«>y>v, which are characteristic, respectively, of terrestial and 
astrophysical conditions, we are able to conclude that while in the former 
case instability always arises as cellular convection, in the latter case it 
will develop either as cellular convection or as over-stability depending 
on the magnitude of the velocity of the magneto-hydrodynamic wave. 
We shall now show how under quite general conditions we can distinguish 
the manner of the onset of instability. 
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As in §§ 8 and 9 we shall again restrict ourselves to the case when both 
surfaces confining the fluid are free. As we have seen, if instability should — 
develop as cellular convection, then the critical Rayleigh number for 
instability is determined by the minimum of the right-hand side of the 
expression (cf. eqn. (84)) 
RO ant Fy (w+ ‘). ree ie 
where for comparison with (141) we have replaced (1-2) in eqn. (84) by y. 
On the other hand, if instability should develop as over-stability then 
the critical Rayleigh number for instability is determined by the 
minimum of the right-hand side of the expression (cf. eqn. 146)) 

Gang) ee S10 Se 2 Q nv(n+v) 
R‘ i a ea poe E Spl ienas (187) 
In eqns. (186) and (187) we have attached superscripts (con) and (0—s) 
to R to emphasize the nature of the instability which follows. 

Using the results of table 1 we can readily draw in the (#,, Q)-plane 
the critical Rayleigh numbers derived on the basis of either of the 
eqns. (186) and (187). Let the curves so derived be called the ‘ convection 
curve’ and the ‘over-stability curve’ respectively. Both are 
monotonically increasing curves. The convection curve starts on the 
#,-axis at the point 


27 ; 
ROM = = at (Q=0),'> lee Cees 
and for @-—co becomes asymptotic to the line 
REM =7?Q (Qo), . . . (189) 
On the other hand, the over-stability curve starts on the R,-axis at 
27 ,8,8.—KV 
(Ce ya ae 7) ees 
R, certo (Q=0), . .°>. (190) 
and for Q— co becomes asymptotic to the line 
(Os) 2 u(n+v) 
Ei Sisear ra (Qc) "Si ee ony 
Since 
Bite Sei ey ae ata sas 9 
Kv(K-+v) (K+) win y) lesa) ean Lye ere ce 


it follows that the over-stability curve always starts above the convection 
curve. Consequently, we need distinguish between only two cases: 

(i) The over-stability curve lies entirely above the convection curve.—From 
eqns. (189) and (191) it follows that this will be the case whenever 


n>. Nr Coe tone ee ee es 


Moreover, it is clear that in this case, when instability arises it will 
always be as cellular convection. This is what happens under normal 
terrestial conditions. : 
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(ii) The over-stability curve intersects the convection curve and for Qc 
lies below it.—According to eqns. (189) and (190) this will be the case 
whenever 


ai ore eer! | 32194) 


Conversely, whenever 7 <x the two curves must intersect. Let the value 
of @ at the point where the two curves intersect be Q*. Then it is clear 
that for @<Q*, instability, when it arises, will be as cellular convection ; 
but for Q@>Q* we shall have over-stability when the critical Rayleigh 
number #,°~ is reached. This is essentially what happens (though 
in a somewhat ‘extreme’ form) under astrophysical conditions. 


§11. ConcLupInG REMARKS 


From the numerical example considered in § 6 (eqns. (113) to (115)) 
it would appear that under normal terrestial conditions (y>«) the 
inhibition of convection by a magnetic field is a very pronounced effect. 
And since convection and turbulence are often postulated in the earth’s 
core (cf. Bullard 1948) it is possible that a study of the problem considered 
in this paper for a spherical geometry may lead to results of some 
significance. The problem of the instability of a fluid sphere heated 
within has recently been investigated by Jeffreys and Bland (1951); 
and this investigation provides an excellent context in which the theory 
developed in this paper can be extended. The consideration of this same 
problem when «> and over-stability is possible may also be of interest 
for astrophysics. However, under astrophysical conditions («>7>v) 
the effects considered in this paper are probably of most direct application 
to sunspot phenomena. In particular, a re-examination of the problem 
of thermal stability of a fluid enclosed in a long circular cylinder which 
has been considered by Hales (1937) in connection with a theory of 
geysers may be instructive. A problem of a different kind relates to the 
nature of turbulence which will finally establish itself in an extended 
atmosphere on account of thermal instability. A beginning in this 
direction has recently been made by the writer (Chandrasekhar 1952). 
An extension of this investigation to the case when a magnetic field is 
present would be of interest. And finally, one might also investigate 
the related problem (cf. Low 1929) of the stability of a liquid between 
two co-axial cylinders rotating at different rates, when the liquid 
considered is an electrical conductor and a magnetic field is present ; as is 
well known this problem in the absence of a field has been investigated 
by Taylor (1923). The author hopes to publish the results of his 
investigation of these various problems suggested by this paper in the 
near future. 

In concluding this paper, I wish to record my indebtedness to 
Miss Donna Elbert for her very valuable assistance in having carried 
out all the numerical calculations involved in the preparation of 


tables 1 to 3. 
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SUMMARY 


The variation in the grain-density of tracks produced in photographic 
plates by electrons movin gat relativistic velocities has been determined. 
It is found that the grain density increases by about 5% between the 
energies 9 and 25, measured in rest mass units, thereafter remaining 
constant at least up to L/mc?~104. 


§1. INTRODUCTION 

PaRTICLES of charge /e/ moving in photographic emulsions with velocities 
close to that of light, produce tracks in which the grain-density is close 
to the minimum value. The precise variation of the grain-density as 
a function of energy is of interest because it may provide a method of 
mass discrimination in certain energy regions, and because it is a problem 
which lends itself to quantitative theoretical calculation. Electron 
tracks are the most suitable for an investigation extending over a wide 
range of values of H/mc?. 

liford G5 emulsions, 400 microns thick were exposed at an altitude 
of 100 000 feet. After suitable processing one 6 in. x6 in. plate was cut 
into four pieces, 3in.x3in. The methods described by Hooper et al. 
(1951), were used in a systematic scan of certain regions, indicated by 
the cross-hatched areas of fig. 1, for pairs of electrons produced by the 
conversion of energetic photons. The close pairs of tracks thus produced 
can be attributed to electrons with great confidence. The energy of the 
pairs of electrons was determined by the multiple-scattering technique 
(Gottstein et al. 1951). In this way, more than 300 electron tracks, 
each at least 2 mm long and with values of H/mc? between 20 and 
1-5 104, were examined. Electron tracks for which 12>EH/mc?>6 were 
obtained by measuring single, highly scattered tracks with an ionization 
close to the minimum value, found in the same region of the plates. 


§2. PLaTE CALIBRATION 
Since the effect to be investigated is a small one, considerable time was 
spent in determining the variation of grain-density with depth and 
position in the emulsion. For this purpose, electron tracks in the region 
E/mc?~500 to 1500 were used, the specific ionization being assumed 
constant in this range. 


* Now at the University of British Columbia, Canada. 
+ Communicated by Professor C. F. Powell, F.R.S, 
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A grain was counted as one unit irrespective of size. This convention, 
useful up to 1:5 the minimum value, has the following advantages :— 

(1) It allows rapid counting, no time being lost in estimating the 
number of grains in aggregates. (2) Subjective errors are small. The 
same individual tracks, recounted weeks later, always agreed with the 
earlier count to within 05%. (3) The variation of the grain-density with 
depth is much smaller than when a convention which * weights * various 
erain-sizes is used. This is presumably associated with the fact that 
grains near the surface of the emulsion are commonly larger than those 
near the glass surface, because of greater physical development. 


Fig. 1 


Cc 


The cross-hatched areas indicate the regions scanned for electron tracks. The 
four 3 in, x3 in. plates were originally part of a single plate. 


The depth variation was investigated for the cross-hatched region of 
plate D in fig. 1. The top 50 and bottom 50 microns of the emulsion 
were not included in this calibration. The remaining 300 microns were 
considered as made up of three 100 micron layers. The grain-density of 
a track in the upper 100 » layer was found to be 0-45% higher, and in the 
bottom layer 0-26% lower than in the centre layer ; the statistical 
deviation for each layer was 1-0% (defined as c=100/1/N, the percentage 
savas deviation, where N is the number of grains counted). With this 
counting convention, there was the igni variati 1 
eS Pie fa, s therefore no significant variation of grain 


° 
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Using this result, a search was made for local variations in the scanned 
areas of plates A, B,C and D. Plate A hada grain-density 1-58°% lower, 
C, 084% lower, and D 2-81% higher than plate B which was taken as 
standard. (For each measurement o was 1:0°%.) The difference of 
439% between plates A and D is apparently significant. It is 
presumably either due to a variation of sensitivity of the poured emulsion, 
to non-uniform development, or to both factors. The variation between 
plates poured and developed at different times is likely to be even 
greater. 


§3. EneRGy DEPENDENCE OF THE GRAIN DENSITY 


The grain-densities of tracks of electrons for the intervals of H/mc? of 
(1) 6-12 ; (2) 20-32 ; (3) 50-150 ; (4) 500-1500 and (5) 5 000-15 000 were 
determined with the help of the preceding calibration. (o was 1-:8%, 
13%, 1%, 1% and 1% respectively.) The results are plotted in fig. 2. 
The grain density is constant to within 0-25° for the energy ranges (2) to 
(5). The decrease of 4:8°% for group (1) is apparently statistically 
significant. It may be noted, that, if the calibration were not applied, 
points (2) to (5) would still be constant to within 1-0% and the decrease 
for group (1) would be 5-5%. The results are in substantial agreement 
with the work of Pickup and Voyvodic (1950), which extends over a small 
range of values of H/mc?. 

The first detailed theoretical calculation of the ionization losses in 
condensed media was made by Fermi (1940). He showed how polarization 
effects modified the Bethe-Bloch theory, making the simplifying 
assumption that the dielectric properties of the medium were 
characterized by only a single dispersion frequency. Wick (1943) and 
Halpern and Hall (1948) put forward a multi-frequency theory. Bohr 
(1948) and Schénberg (1950) have also studied this problem and, in 
particular, have considered in some detail the effect of the Cerenkov 
radiation. All these theoretical treatments are entirely classical, but 
Allen (1951) has investigated the problem from a quantum standpoint. 

Messel and Ritson (1950) pointed out that in the photographic plate, 
it was not the total energy absorbed per unit volume that should be 
considered, but rather only that part which was effective in latent image 
formation (Webb 1948). They therefore suggest that the contribution 
to the absorption of energy due to 5-rays of range greater than the mean 
grain-size, should be neglected. This is done by applying a ‘ cut-off ° 
to the maximum energy equal to 5 or 10 kV. 

There is a lack of experimental data concerning the oscillator strengths 
and excitation energies for the electron transitions of silver bromide, 
particularly in the ultraviolet region. Theoretical calculations on the 
basis of the Hartree method: would lead to only very approximate 
eigenfunctions because of the occurrence of the phenomenon of ‘ over 
stability ’. Consequently, no computations for the ionization loss. in 
AgBr were made. In fig. 2, the curve is a calculation for iron, using 


202 
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a cut-off energy of 5 kV. The corresponding curve for AgBr can be 
expected to differ from this by a few per cent. However, the theory 
appears to be capable of accounting for the main features observed in this 


experiment. 
Fig. 2 


4 
100 { t —e 
95)\- | 

Ls = 2 3 7 


ENERGY (E/Mc* UNITS) 


The observed variation of grain density with H/mc?. The calculated curve is 
for iron, calculated from the theory of Halpern and Hall. 


Fig. 3 
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Probability curves designed to investigate latent image fading. The experi- 
mental points are indicated with maximum reading errors. The grain: 
density is measured per 50-6 microns of track length. 


It would be of interest to investigate the grain density at E/mc2?~2 or 3 
but even with an exposure to a homogeneous electron beam at this 
energy, the experiment is very difficult owing to the large scattering in the 
~ emulsion. A mono-energetic meson beam in this region is not yet 
available. 
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§4. FapiIne oF THE Larent Imac 


The tracks were examined in order to determine if there were 
appreciable fading effects in the G5 emulsion, due to any extreme 
conditions of temperature and humidity during the course of a high- 
altitude exposure. The histogram of fig. 3 shows the distribution in 
grain density normalized with the aid of the previous calibration, of 
57 tracks of electrons whose energies were in the ranges (2) to (5). The 
length of each track, starting from the origin of the electron pair, was 
determined by the convention that it has N=500 grains. Since, for a 
normal distribution, 68% of the tracks should lie within N+-./N and 
N—V\/N, a broader distribution would imply that fading was appreciable. 
To test this point, the experimental data were plotted on probability 
paper (in which a normal distribution becomes a straight line), and the 
result is shown in fig. 3. The arrows indicate the 50th (the mean) and the 
16th and 84th (the limits set by the statistical deviation) percentiles. 
AA is the calculated normal curve, while BB is a better fit for the points. 
Hodgson (1952) has pointed out that because of the finite grain size, the 
statistical fluctuations are somewhat reduced, an effect analogous to the 
counting losses of a Geiger-Miiller counter due to finite resolving time. 
It may therefore be concluded, since there is no significant difference in 
slope between the curves 44 and BB, that any effects due to fading were 
small and without significant influence on the results of this experiment. 
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SUMMARY 


Glauert has recently published photographs showing the distortion 
of multiple-beam interference fringes produced by inserting a thin sheet 
of glass or mica between two interferometer mirrors. It is shown that 
multiple reflections between, the surfaces of the inserted sheet and the 
mirrors necessarily produce the types of distortion observed, and that the 
published photographs can be analysed quantitatively. 


$1. INTRODUCTION 


Tuis paper is concerned with the interpretation of some remarkable 
interference fringes recently published by Glauert (1951). ° The Fizeau 
fringes given by a pair of flats, placed at a small angle to each other, 
are straight, lying parallel to the edge of the wedge. If the flats are 
silvered and the monochromatic illumination is correctly collimated, the 
fringes are formed by multiple reflections and are very sharp... It 
is often proposed to place a (biological) specimen under examination in 
a suitable fluid between silvered optical flats, in the hope that the resultant 
fringes will give some indications as to structure and form. The present 
results show that the greatest care must be taken in the interpretation 
of such interferograms, and they illustrate the ease of production of 
artefacts.’ The simplest results to interpret are those shown in figs. | and 2 
of Glauert’s paper which are reproduced in Plate XXII (a) and(6).¢ In 
Plate XXII(a), a glass cover slip had been inserted between the flats, and 
the resulting fringes were modulated into zigzags, the amplitude of the 
zigzag itself varying along each fringe. The positions of minimum 
amplitude did not correspond on neighbouring fringes, but were 
approximately repeated every third fringe. In Plate XXII (6), a sheet of 
mica had been inserted between the flats. There was a region of violent 
disturbance at the left of the picture, which we do not attempt to 
interpret. Elsewhere the fringes were doubled, and the two components 
could be separated by using polarized light. In one component the 
amplitude of the zigzags was uniformly small, while in the other it was 
uniformly large. Except on two fringes near the middle of the field, 
the zigzags of large amplitude were very unsymmetrical, actually 


overhanging in opposite directions in the fringes on opposite sides of 
these two. 


* Communicated by the Authors. 
+ For plate see end of issue. 


On the Origin of Glauert’s Superposition Fringes 539 


As Glauert suggests, the observed effects can all be explained in terms 
of multiple reflections between the mirrors and the interposed sheet, and 
between the two faces of this sheet. The optical path for a simple passage 
from the upper mirror to the lower is, of course, independent of the 
position of the interposed sheet. In practice a ray cannot make this 
simple passage, but parts of it suffer single or multiple reflections from the 
surfaces of the sheet. The amplitude of the resultant ray depends on the 
phases of those disturbed components relative to the component which 
has made the simple passage, and therefore on the position of the sheet 
between the mirrors. If the sheet is tilted about an axis lying 
perpendicular to the beam of light and to the intersection of the mirrors, 
the effective optical thickness is not constant along one of the original 
fringes, and the fringes observed can no longer be straight lines parallel 
to the intersection of the mirrors. 

For a detailed interpretation, we have to make the reasonable 
assumptions that the glass slip was very nearly plane, but that its thickness 
varied gradually by a few fringes in crossing the field, while the slip of 
mica was of very uniform thickness, but appreciably twisted out of 
a plane. 

In § 2 we discuss the multiple reflection of a single ray between two 
mirrors with a sheet of material of thickness 7’ and refractive index 
u lying between them. In §3 we discuss the fringe pattern formed by 
a plane parallel sheet lying obliquely between the mirrors, and in § 4 we 
consider the additional effects which arise if the faces of the sheet are 
slightly inclined to one another to form a wedge. 

Our interpretation of these fringes suggests that the corresponding 
artefacts will be much less noticeable in the case of a biological specimen 
ina fluid of refractive index close to its own than in the present simple case. 
Firstly, the lines of maximum rate of change of refractive index will in 
general make finite angles with the primary beam of light, and the 
reflected components will therefore not interfere with the direct beam. 
Secondly, the amplitudes of the internal reflections will be small if 
violent changes of refractive index are avoided. 


§2, MuxtreLe REFLECTIONS OF a SINGLE Ray 


We consider an incident ray with electric component A exp (2ka) 
striking the upper silvered surface (x0), suffering multiple reflections, 
and finally emerging as a reflected ray B exp (—ikx) and a transmitted 
ray W exp (ikx). Between the upper silvered surface and the upper 
surface of the glass plate, the reflected rays combine to form resultants 
C exp (ikx) propagating downwards, and Dexp(—vtkx) propagating 
upwards, and similarly in the other two regions (fig. 1). The refractive 
index of the glass is ., and the distances X, Y and 7 are indicated in the 
figure. 

Let t and r be the coefficients of transmission and reflection at a glass-air 
interface, defined so that the transmitted flux of energy is |t|? times the 
incident flux, while the phases of the electric vector in the transmitted 
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and reflected waves are given correctly by ¢ andr. Then for a surface 
at 7=0. 


t=2u1?/(u-1), | 
r=+(u—1)/(u+1), Petite (1) 
the sign depending on the direction of the incident ray, 
and 72+ f=], 
Fig. 1 


|eeie f Retkx 


[ce ZK 2 { pEetkx 


| Weeks 


Multiple reflections of a single ray. The incident ray A gives rise to a reflected 
ray B and a transmitted ray W. 
The corresponding quantities for the silvered mirrors, ¢/ and 7’ are in 
general complex, and |?’ |?+-|r’|?<1. We can, without loss of generality, 
take the mirrors to be thin silvered sheets with no glass backing. 
The continuity of E at each interface leads to the equations 


B=r'A+4t'D, 
C=t'A-+r’'D, 
De7 2-8 VY — CoM —2-¥) 3. Lite Wha = T= Fy 
Fett\< -T-¥) 1-114 ea— TY) + pGle~#X-T-¥) 
: j > 
J Bos ey Mee a “ x é 9 
Ge ux =e )— pp Peth(X—-I 4 y12¢Ne— MH X—I ) ( ) 


Me X~Y) 1124 PeieX—¥)_pNe-i(X-¥), 
Ne—*<—r' Me*<, 


and W=t'M. 
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Here r stands for +(—1)/(u-+1). 
Solving these equations for the ratio W/A, we find 


| W/A |=t?2/n, ean e 29- Os ees (3) 


where 7 is a rather complicated expression. We know that as the thickness 
of silvering on the mirrors is increased, so that |7’}+1 and |#’|>0, the 
interference pattern becomes fainter, while the bright fringes become 
narrower. In the limit |7’|+1, the pattern disappears completely, but 
the positions of these fringes are given by the condition 70, which can 
now occur for real values of X and Y. By adding constants to X and 
Y, we may take r’= +1 for both mirrors. The condition 70 then 
becomes ° 

(1p re~2#¥) pe 2iukl (y» 1 9 2ik¥) 

=O ee WEE (pag 2 Ey. 4(] bre 2*Y)T. May (4) 


The ray considered contributes to a bright fringe if this condition is 
satisfied, and if X, Y and 7’ are functions of position, the equation of a 
fringe is 

X—T=p'(Y)+(2n+1)zx/2k 
(l—r?) cot pkT+2rsnm2kY;-, . . . . (8) 


where ~ tan kp’= (1-++-r?)-+ 2r sin 2kY 


' which may also be written 


X+(p—1)T=p(¥)+na/k 
r sin 2ukT'+ 2r sin pkT cos (2k Y +pkT’) | ey CAE) 


where tan ke=— 773 cos 3b T}2r sin pT sin (EY +pkT) 


Here n is any integer. 


§3. Frrvees Formep By INCLINED MIRRORS WITH A SHEET OF UNIFORM 
THICKNESS [INTERPOSED 


We now suppose that the various surfaces in fig. 1 are not exactly 
parallel, although so nearly parallel that multiply reflected beams will 
combine to give interference. In the plane of the lower mirror we choose 
axes wu, parallel to the intersection of the mirrors, and v, perpendicular 
to this (fig. 2). Then X is a linear function of v 


X=6v. Rapparees ier cany 9-7 (7) 
If we take »—1, r=0, the fringes are given by (6) as 

Tn] Che sa 

p=0. 


Whatever the value of r, the numerator in the expression for tan kp is 
bounded, while the denominator obviously cannot vanish if 1—r?—2r>0 
ie. r<1/2—1,u<1/2+1. (The denominator is, in fact, always positive 
if r<1/4/2, w<3+2,/2.) For any practical value of y, tan kp remains 
bounded as ,. increases from 0 to its actual value. If we number the 
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fringes in the absence of a glass plate by the integer n, we may establish 
a 1:1 correspondence between the fringes and those in the presence of 
a glass plate by adding to (6) the condition 


—dn<kp <hr. 


Fig. 2 


fas 


(0) 
Sections showing the glass plate between the two mirrors. (a) section normal to 
the intersection of the mirrors, (b) section parallel to the intersection of the 
mirrors and to the rays of light. 


Suppose now that the glass plate is not parallel to the intersection of 
the mirrors, so that Y is a function of wu. It follows from (5) or (6) that 
the fringes are no longer straight lines, v—constant, but are modulated 
in the v direction. For a given thickness 7’, the amplitude can be 
determined from (5). If we put tan kY=—Z, there results a quadratic 
equation in Z, the roots of which coincide if 
1+r? 
1—r? 


2r 
tan kp’ = 


cot wkT' + 


cosec pkT’. 4 Ee) 


l—r? 
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If we call the corresponding values p,’ and p,’, we have 
tan k(p1'—p2')=2¢/(1—¢?), 


o=2r sin pkT/(1—r?). 


where 


Since the extreme values of p’ must coincide if r—0, the appropriate 
solution is 
k| py'— ps’ |=2 tan-1|2r sin pkT/(1—r?)|. . . . (10) 


This gives the total amplitude of displacement of a fringe |p,’— Pp, | on 
the scale for which the distance between neighbouring fringes is z/k. (The 
same formula (10) may be obtained from (6) : there is then no ambiguity 
in the choice of solution.) 

In the case of the glass cover slip, the thickness varies appreciably 
in the field, and the amplitude varies from zero to its maximum value 
of (2/7) tan~1|27/(1—r?)| times the spacing of fringes. By measuring 
the maximum amplitude and the spacing, we, find 2r/(1—r?)— tan 26-2°, 
giving 41-61, a reasonable value. In the case of mica, the thickness is 
apparently constant to a small fraction of a fringe. This thickness is such 
that for one polarization sin ukT is practically zero. If we assume that 
for the other polarization sin .kT is unity (so that the mica formed a 
quarter-wave plate), we find »=1-57, while the principal refractive 
indices of mica lie between 1-56 and 1-60. 

We now consider in more detail the effect of inserting a slightly tilted 
plate of uniform thickness. In addition to eqn. (7) for X, we have 


Deed OU EDO, itl he wien een? Sete Ld) 


_ where « and f are small angles comparable with 6. The equation of the 
fringes in the (u, v) plane is given by (6) as 


dv-+(u—1)T'=p(Yotau+tfv)+na/k —4nr<kp<}n. . . (12) 


It follows from (6) and (12) that a given fringe undergoes one complete 
cycle of modulation when kaw increases by 7, while we pass from one 
fringe to the next when k6v increases by 7. The observed ratio of the 
distance between fringes to the wavelength of the modulation gives 
a/@ directly : for the glass cover slip we find «/9=1-73. For the mica, the 
passage from fringes of short wavelength at the top of the figure to those 
of long wavelength at the bottom indicates that the sheet was warped, 
dipping (or rising) steeply to the right near the top, and less steeply near 
the bottom. 

The form of the modulation, given by (6), is noticeably unsymmetrical. 
It is sketched in fig. 3 (a) for the case sin wkT'=1, w=1-61. 

The effect of the term fv is to shear the whole pattern through an 
angle —/«, the lines v=constant remaining parallel to themselves 
(‘simple shear’). As is shown in fig. 4, this shear may turn the original 
modulations into overhanging teeth. This is actually observed in 
Plate XXII (b). The fringes at the top of the figure overhang to the left, 
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those in the middle are symmetrical, and those at the bottom lean to the 
right. We deduce that the mica was bent about an axis roughly parallel 
to the fringes. The fringes observed with glass do not overhang, and are, 
in fact, more nearly symmetrical than those of fig. 3 (a). We deduce 
that they have been sheared through an angle of about tan~! (0-4) when 
their form should be as in fig. 3 (6). It follows that f/@=-0-7. 


Fig. 3 
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Calculated shape of the fringes of maximum amplitude, with refractive index 
=1-61. (a) undistorted, (b) sheared through an angle of tan-1 0-4. 


§ 4. Krrecr or A GRADUAL VARIATION IN THE THICKNESS 
OF THE SHEET 


We shall now show that three other effects observed in the fringes 
formed with a glass cover slip can be explained on the assumption that 
the thickness of the glass 7’ was not constant, but varied linearly in the 
field of the picture according to the law 


[as Te, leis Oneien Mate. 4h ac Lo 
These effects are (i) that the mean spacing of the fringes where the glass 
is interposed is not equal to the mean spacing of the undistorted fringes, 
(ii) that the amplitude of modulation is not constant along any one fringe, 
and (iii) that points of minimum amplitude on neighbouring fringes do 
not correspond, while those on every third fringe usually correspond. 
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We again assume that 7’ varies so slowly that the glass in any 
neighbourhood may be considered as of uniform thickness. The equation 
of the fringes is then 


du-+(u—1)(To+yu+dsv)=ptnnfk . . . . (14) 
Fig. 4 
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(b) 
Effect of the shear produced by the angle of tilt « (schematic). (a) fringes 
before shearing, (b) fringes after shearing. 
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and the condition —47<kp<47 ensures that the term in p does not affect 
the mean direction or spacing of the fringes. 

The mean direction of the fringes is such that [@+(u~—1)d]vu+(u—1)yu 
is constant, while the direction of the undistorted fringes is v—constant. 
Thus the mean direction of the fringes turns on meeting the glass through 
an angle tan—!{(~—1)y/[0@+(~—1)6]}. This cannot be observed in the 
picture. The mean spacing of the fringes in the v direction is 
m/k{@-+-(u—1)8]. In the absence of the glass the mean spacing is 7/k@. 
The ratio of these spacings is 1+(u—1)6/6 : the observed ratio is 8/7, and 
we have determined »=1-61, so that 6/@=0-2. 

The variation of amplitude along a fringe is given by (10). The 
modulation vanishes whenever pkT’=nr. If the points (w,,v,) and 
(Us, Vy) lie on the same fringe, we have from (14) that 


A(¥y — V2) + (H—1)[y(Uy — Ug) +8(Y, — V2) ]=0, 
while if they are neighbouring points of zero modulation we have 
wk y(u— Us) +8(v) V2) =n. 
Eliminating v, we find 
(uy —Ug) = + (u— 18/6] /ehy. 
The ratio of this to the mean spacing of fringes in the v direction is thus 
| [1+ (u—1)8/OPO/py. 


The observed ratio is about 11. Since [1-++(—1)3/0]=8/7, we have 
y/@=0-076. The angle through which the fringes should be bent on 
meeting the glass is thus tan-! 0-:04<24°, and it is reasonable that this 
is not observed. 

Finally we check this portion of the analysis by means of the 
observation that the positions of zero amplitude on roughly every third 
fringe correspond. For we go from one zero to another by moving in 
the v direction a distance z/uk5, which is (6/ud)[1+(~—1)8/6] fringes. 
Using the values of 6/9 and w already determined, this is equal to 
3-03 fringes, agreeing satisfactorily. 
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XLIX. An Unusual Birefringence Interference Pattern 


By 8. Tontansky and F. 8. A. Suuran 
Royal Holloway College, Egham* 


[Received February 11, 1952] 


In a recent note to Nature, Glauert (1951, Nature, Lond., 168, 861-2) 
has described, but did not explain, some unusual localized interference 
fringes obtained by placing a thin sheet of mica between two suitably 
silvered optical flats. On viewing with monochromatic radiation a 
highly complex system of localized multiple beam fringes form, with 
modulations arising through the wrinkling of the mica which, in fact, 
leads to a triple interference film, two air gaps and one mica gap. 
A striking feature of the fringes described by Glauert is the anomaly 
arising from the birefringence of the mica. The fringes are double, 
the two components being mutually perpendicularly plane polarized, 
but they do not follow identical modulation patterns, being variously 
straight or zig-zag according to location, even at times when very close 
to each other. The modulation effect is clearly determined both by 
the geometry and the refractive index. 

We wish to describe here another unusual interference phenomenon 
which is closely allied to and which throws some light on the behaviour 
of the fringes described by Glauert. A drop of molten para-toluidine 
was placed between two silvered optical flats and allowed to crystallize 
in situ. The toluidine contacted the surfaces and therefore there were 
no effective air films. On examination with mono-chromatic light (green 
mercury) the system showed a distribution of multiple beam Fizeau fringes 
which were double owing to the crystal birefringence. One small region 
which is illustrated in fig. 1 (a) (Plate X XITI}) exhibited an anomaly which 
we have not seen elsewhere in any of the many hundreds of interferograms 
recorded in this laboratory. The fringe pattern consisted of a network of 
two systems mutually plane polarized perpendicularly, these two systems 
crossing each other almost at right angles. Fig. 1 (a) shows the fringe 
pattern whilst 1 (b) and 1 (c) show the two patterns obtained respectively 
by the introduction of polaroids, set, in the two cases, at right angles to 
each other. This system of two independent sets of fringes can only be 
interpreted one way. The crystal is acting as a simple wedge to each 
direction of polarization, but the two ‘ wedges ’ are almost at right angles. 
One is confronted with a position in which there is a geometrical wedge 
in one direction and an optical wedge in another direction. The 
geometrical wedge forms because the surfaces of the flat are not parallel. 
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The crystal axes must be so oriented in this particular form of 
crystallite such that with respect to the two different polarizations the 
effective wedges are roughly at right angles to each other. 

This curious fringe pattern appears to be related in form to those 
described by Glauert, but much simplified through the absence of the air 


films, and the wrinkling. 
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L. On the Significance of the Limit Load Theorems for an Elastic-plastic 
Body* 


By E. H. Lest 
Graduate Division of Applied Mathematics, Brown University, U.S.A.t+ 


[Received January 15, 1952] 


ABSTRACT 


The limit load theorems of Drucker, Greenberg, and Prager (1951) for an 
elastic-plastic material flowing under constant surface tractions are 
contrasted with the theory of the yield of a plastic-rigid body (Hill 1951). 
Plastic-rigid theory is interpreted, as it must be for application to actual 
materials, as a model providing an approximation which can reasonably be 
expected to be satisfactory when plastic strains large in comparison with 
elastic strains have occurred. It is shown by the discussion of particular 
examples that the limit load for an elastic-plastic body is in general reached 
or closely approached while strains throughout are of the order of elastic 
strains, and at this stage the limit load theorems can be applied. Examples 
are cited to show that if plastic-rigid analysis is adopted, and strains large 
compared with elastic strains must be attained, the influence of boundary 
motion may invalidate application of the yield limit theory in circumstances 

_when the limit load theorems for an elastic-plastic material provide 
significant results. 

It is shown that the uniqueness and variational theorems applicable to a 
plastic-rigid material at yield, apply to an elastic-plastic material at the 
limit load, thus permitting stress distributions in an elastic-plastic body 
to be determined on the basis of plastic-rigid theory. Examples are cited. 


$1. INTRODUCTION 


WE are concerned with bodies of elastic ideally-plastic material, having a 
stress-strain relationship of the type shown in fig. 1. Quasi-static motion 
is considered for which inertia forces can be neglected. When such bodies 
are subjected to external tractions, as a rule they are initially stressed 
elastically throughout. With increasing tractions plastic flow will 
commence, but there will in general be a stage of contained plastic flow, in 
which the regions flowing plastically are constrained by surrounding 
elastic material so that strains throughout are of the order of elastic 
strains. With increasing loads, the plastic regions will in general spread 
until such a constraint ceases to apply, and continued flow takes place with 


* The results presented in this paper were obtained in the course of research 
sponsored by the Office of Naval Research. 
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little or no increase in load. This situation is termed uncontained flow, 
and the constant load at which continued flow can occur is termed the 
limit load. Recently Drucker, Greenberg and Prager (1951) have developed 
theorems prescribing whether or not such continued flow will take place 
under specific boundary conditions (see also Prager and Hodge 1951, 
chap. 7). When particular situations are considered, these theorems can 
be used to determine upper and lower bounds to the limit loads. The 
power of these theorems lies in the circumstance that they provide 
information about the uncontained flow, without the necessity of analysing 
the previous contained flow which is in general beyond the scope of 
presently available methods of analysis. ’ 
Recently Hill (1951) has discussed the state of stress at the yield point 
of a rigid-plastic body, that is, from the standpoint of application, a body 
in which elastic strains are neglected so that the elastic line in fig. 1 is 
considered to lie on the o axis, deviation from it occurring only at the 
yield point stress Y. Such a hypothesis can only be expected to provide 
a satisfactory approximation for an actual elastic-plastic material if the 


Fig. 1 


PLASTIC FLOW. 


=-ELASTIC DEFORMATION 


Stress-strain relation. 


plastic strains occurring are large compared with elastic strains, so that 
the latter can be neglected compared with the former. For a rigid-plastic 
material, no deformation occurs until uncontained flow takes place, and 
Hill presents theorems concerning this initial yield which are in many 
respects similar in form to the earlier theorems of Drucker, Greenberg and 
Prager for an elastic-plastic body. Hill’s theorems are obtained by a 
generalization and reinterpretation of his work, Hill (1948, 1950), on the 
principle of maximum work and its complement, and he states that the 
theorems of Drucker, Greenberg and Prager are special cases of this 
generalized principle. To show that this is not the case is one of the objects 
of the present paper, and it is demonstrated in the next section that in 
certain cases it is essential to consider the deeper elastic-plastic analysis in 
order to be able to make useful application of such theorems in practice. 

The reason for this difficulty is that in the application of both groups of 
theorems, the initial unstressed surface of the body is used to specify 
boundary conditions. For a plastic-rigid body at yield, previous strains 
being of the order of elastic strains are logically neglected, and the unde- 
formed surface is therefore used, For the elastic-plastic body the 
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theorems are expressed in terms of the current surface geometry when 
uncontained flow is impending, but since in general the contained flow has 
not been analysed, the deformation of the boundaries is not known, and so 
the initial undeformed boundaries must be used as an approximation to the 
current boundary geometry. For many problems involving strains of the 
order of elastic strains the motion of the boundaries has a negligible 
influence, in fact the usual linear theory of elasticity is based on this hypo- 
thesis. | However, as the slenderness of the geometrical form increases, 
boundary motion has an increasingly important influence, and for example 
in the case of elasticity we reach the group of elastic instability problems 
for which boundary motion dominates the solution. The influence of 
boundary motion clearly becomes more important as the strains increase, 
so that application of a rigid-plastic analysis, which can only be expected 
to give useful results for an actual body when strains occur which are 
large compared with elastic strains, may imply a marked influence of 
boundary motion. An analysis neglecting this will then be unsatisfactory. 
In many cases the elastic-plastic analysis leads to attainment or close 
approach to the limit load while the strains are still of elastic order, 
say twice the elastic limit strains, so that the limit load theorems can 
be used to give significant results before boundary motion has an appre- 
ciable influence. Thus the limit load for an elastic-plastic body may be 
significant, whereas the practical application of the corresponding plastic- 
rigid analysis which requires strains large compared with elastic strains 
may be vitiated by the influence of boundary motion. Examples are 
given in the next section to illustrate this situation quantitatively. 

Before turning to the discussion of specific examples let us consider the 
first theorem of Drucker, Greenberg and Prager (1951) which contains the 
basic difference between their analysis and that of Hill. It states that if 
changes in geometry which occur during uncontained flow under constant 
surface tractions are neglected, then the stresses remain constant, and so 
only increments in the plastic components of strain occur. ‘Thus the elastic 
strain increment components are zero, and it is for this reason that the 
resulting theorems for an elastic-plastic body bear a close resemblance to 
corresponding theorems for a rigid-plastic body in which from the standpoint 
of practical application the elastic components of strain are simply 
neglected in comparison with the plastic components. 


§2, BEHAVIOUR IN UNCONTAINED FLOW 


Without a detailed knowledge of each particular case the term uncon- 
tained plastic flow, that is flow not constrained by surrounding elastic 
material, is difficult to define precisely. From the study of particular 
examples it is found that sometimes the spread of the plastic region reaches 
a stage at which continued plastic flow occurs under constant surface 
tractions, if geometry changes are neglected, while the strains throughout 
are of the order of elastic strains. In other cases, although there is no 
direct constraint due to surrounding elastic material, for continued flow the 
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surface tractions must increase towards a limiting value, and in general 
almost attain this value while strains throughout are still of elastic order of 
magnitude. Examples of such behaviour have appeared in the literature 
(e.g. Hill, Lee and Tupper 1947), and such situations have recently 
been emphasized by Hill (1951). The conditions under which this increase 
occurs may be thought of as a type of constraint due to the elastic properties 
of the material, for in the notation of the theory covering a particular 
example given in the Appendix, it represents a lag in the orientation of the 
stress vectors behind the strain increment vectors which prevents exact 
‘free’ flow, that is flow under constant stress. Such a lag is due to the 
elastic components of strain. Thus generalizing our concept of elastic 
constraint in this way, limit flow occurs under constant or almost constant 
surface tractions, and to this extent the limit load theorems are exact or 
nearly so. From the study of particular examples it is clear that in general 
for reasonably steady boundary conditions this lack of preciseness in the 
attainment of the limit load concerns differences of very small magnitude 
which are of no importance in the practical application of the theorems. 


Fig. 2 
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Stress distribution in twisted tube. 


Since only a few elastic-plastic problems are amenable to complete 
analysis, we can use only simple problems for illustration. However, a 
quantitative study of these can be used to indicate the type of behaviour 
to be expected in more complex cases. As an example of flow at constant 
limit load, consider the torsion of a thick-walled tube of internal radius 
a, external radius 6. When a torque below the elastic limit is applied 
the shear stress is proportional to the radius and the torque increases 
linearly with the twist. When the shear stress at the outside reaches the 
limit value for plastic flow, k, plastic flow commences there, and with 
increasing torque the plastic region spreads towards the centre, the 
elastic-plastic boundary being a cylindrical surface concentric with the 
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tube surfaces. Figure 2 shows the distribution of the shearing stress 
when the plastic elastic boundary is of radius p. With increasing torque 
the plastic region spreads, and the stress at each point in the plastic region 
remains constant. The plastic-elastic boundary reaches the interior 
surface at a twist b/a times the twist at the elastic limit, and thereafter the 
tube is fully plastic. While part of the tube is still elastic, increasing 
torque is required to produce increasing twist, but when the tube becomes 
fully plastic no further increase in torque is needed. We then have limit 
load conditions and stress rates are zero throughout the tube. For this 
particular case the deformation in the form of twist does not alter the shape 
of the tube, so that we have a special case in which geometry changes do not 
influence the limit load. The limit load theorems of Drucker, Greenberg 
and Prager will therefore apply exactly as soon as the tube becomes fully 
plastic, which for moderate ratios of b/a will be when the maximum strain 
is still of elastic order. 

As an example for which, as far as elastic constraint is concerned, the 
limit load theorems apply exactly for strains of elastic order of magnitude, 
but which involves a boundary geometry influence, let us consider the 
expansion of a spherical shell under internal pressure. This problem has 
been fully discussed by Hill (1950b, chap. V). Lamé’s well known 
solution applies for internal pressure below the elastic limit. Plastic 
flow commences on the inside and at any stage the plastic-elastic boundary 
is a spherical surface concentric with the inner and outer surfaces. With 
increasing internal pressure a fully plastic state is reached when the 
plastic-elastic boundary reaches the outside surface, and thereafter 
continued flow occurs under constant internal pressure if boundary motions 
are neglected. Thus the conditions of the limit load theorems apply as 
soon as the spherical shell becomes fully plastic. 

In fact, because of symmetry, the stress is given by the two components, 
radial and circumferential, which are statically determined by the equation 
of equilibrium and the yield condition when the boundaries are prescribed. 
Taking into account the influence of boundary changes, the internal and 
external radii increase due to the internal pressure, the thickness of the 

shell is reduced and continued fully plastic flow is deduced under decreasing 
- internal pressure. For moderate ratios of internal to external radius the 
above limit load conditions pertain while the strains throughout are of the 
order of elastic strains. In practice of course, once the fully plastic condi- 
tion arises with a non-hardening material, flow under falling pressure 
would be unstable and local thinning would be expected to take place, but 
the condition immediately following the development of full plasticity will 
be correctly analysed by the elastic-plastic limit load theorems. Adopting 
a strictly plastic-rigid analysis, one cannot be assured that deductions have 
any validity for an actual material until strains large compared with elastic 
strains occur, and at such a stage it is to be expected that local instability 
will have appreciably modified the solution. Thus whereas the elastic 
plastic approach provides definite information, nothing can be deduced 
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with certainty from the plastic-rigid theory. It is this opportunity to 
apply the theorems of Drucker, Greenberg and Prager much earlier in 
the development of plastic flow, than the plastic-rigid theorems of Hill, 
which permits the much wider use of the former, in cases where the latter 
on their own cannot be expected to give applicable results. However, it is 
clear from the elastic-plastic analysis, that Hill’s theorems, which are. 
formally identical with those for the elastic-plastic body at the limit load, 
do in fact apply much earlier in the development of plastic flow than can 
reasonably be expected from a plastic-rigid analysis. This arises since 
elastic strain increments become negligible compared with plastic strain 
increments long before the total elastic strain can be neglected compared 
with the total plastic strain. Such an application of plastic-rigid theory 
has been made by Drucker, Greenberg, Lee and Prager (1951) to obtain 
information concerning boundary value problems for an elastic-plastic 
body. It is shown there that by substituting stress and strain distributions 
based on plastic-rigid analysis, the limit load theorems for an elastic- 
plastic body determine the exact limit load for this case. Thus plastic- 
rigid theory is shown to have a specific significance in elastic-plastic 
analysis, quite apart from the approximation reasonably to be expected — 
for an elastic-plastic solution involving large strains. 

Let us now consider an example in which the limit load is approached 
but not reached exactly, even neglecting the effects of boundary motion. 
We consider a tensile test in plane strain in which uniform tension o,, 
produces uniform stress in a rectangular block of material, motion is 
prevented in the z direction so that e, is zero, and the faces normal to the 
y direction are free so that o, is zero. The resulting stress-strain relation 
has been given by Hill (1950 b, chap. IV) but in order to facilitate the 
assessment of the approach to uncontained flow, a solution is given in the 
Appendix based on the plane stress and strain diagram introduced by 
Reuss (1930) and discussed and utilized by Hill, Lee and Tupper (1947). 
_ This solution shows that with increasing tension o,, the block is initially 

entirely elastic before the elastic limit is reached at which the entire body 
becomes plastic. However, because of the difference in Poisson’s ratio for 
elastic and plastic. strains, the stress components at this stage are not 
compatible with free plastic flow in plane strain, and as this state 
is approached g,, rises, and with increasing strain it rapidly approaches an 
asymptotic value, the fimit stress. This is shown in fig. 3 where the 
percentage disparity of o,, below the limit stress is plotted against the strain 
in units of the elastic limit strain. It is seen that o,, is within 0-1°%, of the 
limit stress when the strain is only 24 times the elastic limit strain. Taking 
Ex to represent the time, the rate of increase of surface traction becomes 
insignificant at this stage, and the elastic-plastic limit load theorems apply 
to a high degree of accuracy according to the development leading to 
Theorem I of Drucker, Greenberg and Prager. The stresses remain 


effectively constant throughout, and only plastic increments of strain need 
be considered. 
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For this example we can obtain the influence of boundary motion on the 
_tensile force required for continued extension. The cross-sectional area 
of the block is reduced due to the lateral contraction (—e,), and this effect 
is also shown in fig. 3 as a percentage effect reducing the limit load below 
that based on the neglect of this influence. It is seen that this boundary 
influence increases with the strain and soon dominates the variation in 
surface traction as the asymptotic limit stress is approached. For this 
example the boundary motion effect only attains a magnitude of about 4°, 
for strains of the order 5 times the elastic limit strain, and so is itself 
quantitatively unimportant. However the possibility of an unstable 
necking contraction has not been considered which for long specimens 
would be likely to produce a much larger influence at such average strain 
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magnitudes. The solution for compression in plane strain is essentially 
the same with reversed signs for the stress and strain components, so that 
the increase of the cross-section produces an increase in the limit load. 
However a complete analysis would have to include the influence of 
lateral instability, which would involve a larger boundary motion influence. 
Thus a study of this example indicates that the asymptotic approach to the 
limit stress may occur so rapidly that it is effectively reached for strains 
of elastic order, and the boundary motion may produce a larger influence 
on the limit load at this stage. If larger strains must be permitted before 
the analysis is valid, as in the case of approximating an actual body by 
a plastic-rigid model, the neglect of boundary motion may involve a 
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serious error. For most problems, in which plastic flow occurs locally in 
regions adjacent to elastic regions in contrast to the uniform stress situation | 
of this particular example, the redistribution of stress during contained 
flow is likely to anticipate the asymptotic approach to the limit stress, and 
so make it more rapid, while the concentration of boundary deformation in 
the active plastic regions is likely to produce a more marked influence of 
boundary motion. The latter is of course strongly emphasized in the case 
of slender bodies, and bodies of open section. There will, of course, be 
cases where the asymptotic approach to the limit load will be slower than 
for the particular example considered. Thus in cases where stress discon- 
tinuities occur in plastic-rigid analysis, which must be considered as the 
limit of narrow elastic regions (Lee 1950), slower approach to the limit 
is to be expected. The torsion example considered above for a solid shaft 
(a=0) would be such an example. However, in many cases effectively 
constant surface tractions occur at strains of elastic order while the 
influence of boundary motion is still small. The elastic-plastic limit 
load theorems will then provide good approximations to the conditions at 
which deflections begin to increase markedly without appreciable change of 
load, so that there is a well defined knee in the load deflection relation. 
In many cases thus satisfactorily analysed for an elastic-plastic material, 
the use of a rigid-plastic theory, which demands strains large compared 
with elastic strains before its predictions can reasonably be expected to 
comprise an approximation for an actual material, may not be valid due to 
the large boundary motion influence at such strains. It is for such cases 
that the application of the elastic-plastic theory is essential. 


§3. UNIQUENESS THEOREMS 


If we make use of the development leading to Theorem I of Drucker, 
Greenberg and Prager which shows that at the limit load for uncontained 
flow under constant surface tractions only increments in the plastic 
component of strain occur, the uniqueness and extremum principles proved 
by Hill (1951) for a plastic-rigid body can be shown to apply to an elastic- 
ideally plastic body at the limit load. The proofs carry through exactly as 
in the plastic-rigid case. | 

Thus if (o;, w;) and (o*;,, w,*) are stress and velocity fields satisfying 
boundary conditions prescribing velocities over part of the surface and 
tractions over the rest, and each satisfying equilibrium and _ plasticity 
equations in its own regions of flow, and not violating the yield condition in 
regions of flow of the other solution, then o,; and o* ij are identical through- 
out regions of non-zero strain rate of either w; or u*,. 

Also the principle of maximum plastic rate of work over parts of 
the boundary where velocities are prescribed and its complement apply for 
an elastic-plastic material at the limit load. This compares the rate of 
work for an actual solution, and that for the actual velocities and any 
equilibrium stress distribution o*,, which does not violate the yield condi- 
tion. In the case when tractions are prescribed over the whole surface, 
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this implies that o,, for the actual solution is equal to o* 
deformation is occurring. 

These principles for an elastic-plastic material at the limit load are much 
stronger than can be proved for an elastic-plastic material in the region of 
Errnined plastic flow. 

The fact that at the limit load only plastic components of strain 
increment occur, enables plastic-rigid solutions to determine the stress 
distribution in the regions flowing plastically of an elastic-plastic body. 
The several examples which were presented (Drucker, Greenberg, Lee and 
Prager 1951) to demonstrate exact determination of limit loads, therefore 
also present the exact stress distribution in the plastic regions at the 
limit load for elastic-plastic materials. This conclusion provides a 
fuller means of assessment of the development of numerical solutions for 
contained flow. 

In the proof of the theorems for a rigid-plastic material, the orthogonality 
between strain and stress increments is not needed so that the results go 
through formally for a work-hardening material. This is in contrast to 
the elastic-plastic analysis. However, since for application to an actual 
material, the rigid-plastic theory can only be expected to be satisfactory 
when strains large compared with elastic strains have occurred, the work- 
hardening due to these strains is neglected in the theory of the yield limit, 
and stress magnitudes will be uncertain to this extent. Thus the practical 
extent to which a work-hardening material is analysed is subject to this 
limitation, and within the uncertainty of such stress magnitudes, the 
elastic-plastic theory also applies, since prior work-hardening can be 
included as a variable yield limit. However, such an analysis cannot be 
considered a satisfactory treatment of the deformation of a work-hardening 
material. 


ij wherever 
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Tue TENSION TEST IN PLANE STRAIN 


We consider a block of elastic-plastic material stressed homogeneously 
in plane strain by the tensile stress o,, with «,=0 and o,=0. Discussion 
of the stress strain response is facilitated by using the plane diagram 
obtained by projecting the three dimensional stress and strain spaces 
[plotted with the same coordinates (o,, o,, a,) and (2Ge,, 2Ge,, 2Ge,) where 
G is the elastic shear modulus] on the plane o,+0,+o,=0. Previously 
this has been considered as a plot of the stress and strain deviators along 
axes inclined at 120° to each other, but if we note that equal vectors along 
each axis produce a zero resultant, we can plot the stress and strain 
components directly and simplify the use of the diagram. Figure 4 
shows this stress, strain plane, in which the Mises yield criterion for an 
ideally plastic material takes the form of a circle of radius Y, the yield 
stress in tension. For purely elastic strains, the stress and strain points 
coincide and lie within the yield circle. When plastic flow is occurring the 
stress point must be on the circle, while the strain point can move outside it. 
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It is convenient to specify the stress point during plastic flow by the angle 6 
which the stress vector makes with the horizontal, and the strain point by 
coordinates (£, 7). Then for the particular problem under consideration 


we have : 


Gy=2Y]/r/(8) cos 0, o,= Y/+/(3) cos 6 Y sin @;.” “alice 
E=4/(3) G(e,—€,), N=GlEen +). rea 

The compressibility condition becomes 
(3K/G) n= Y[v/(3) cos 6—sin 0] Ss A eee 


The plane diagram for stress and strain. 


and the Prandtl-Reuss relation connecting the total strain increment 
associated with an increment of stress : 
tan 0=(dyn— Y cos @ dé)/((dé+Y sinédé), . . . . . (4) 


(3) and (4) combine to give a directly integrable differential equation in 
€ and 6, which in combination with the initial conditions based on the 
initial elastic strain determines : 


E= Y {cos 0)— | log tan (0/2) + (G/3K)[1/(3) sin 0+ cos 6+-log tan (6/2)] ee 
(5) 


4) is the direction of the stress vector during purely elastic deformation, 
which by use of (3) and the coincidence of the stress and strain points is 
given by 


tan O95=[4/(3)/{1+-(3K/G)}]. » 2. . e686) 
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In terms of (€, 7) the strain components are given by 


2Ge,= M+lEV(3)]},  G@ey=(n—[é/V(3)]}- - . (7) 


It is clear from fig. 4 that increasing strain <,, corresponds to increasing €, 
and according to (5) this is associated with decreasing 6. In fact as 6 
approaches zero, € approaches infinity, so that as the extension proceeds 6 


Fig. 5 


Ex 


See 


() 


decreases towards zero. In fig. 56 is plotted against the extension «,, 
expressed in units of the yield point extension, and this curve emphasizes 
the rapid decrease of 6 towards zero. Thus from (1), o, approaches Y/4/(3) 
the disparity being of order 0, and o,, approaches 2Y/\/(3) with a disparity 
of order 62. This increase of c,, and associated rapid approach to a limit is 
the asymptotic approach to the limit stress discussed in § 2 of the paper. 
The percentage disparity below the limit stress is plotted in fig. 3 against 
the extension. It is seen that the stress is within 0-1% of the limit stress 
when the strain is only 2} times the elastic limit strain. 
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In order to obtain the variation of total load the increase in area normal 
to o,, must be considered by computing ¢, from (7). This can be expressed 
as a percentage reduction of the tensile load below that based on the initial 
area. For this simple case, this represents the influence of the motion of 
the boundaries on the limit load, and is plotted in fig. 3 to provide a com- 
parison with the influence of the asymptotic approach to the limit stress. 
It will be seen that for strains greater than three times the elastic limit 
strain, the reduction of area influence, although small, completely 
dominates the influence of the asymptotic approach to the limit stress. 
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Comments on Dr. Lee’s Paper* 


By R. Hm 
Department of Theoretical Mechanics, University of Bristol 


Dr. Lee claims that the so-called plastic-elastic limit theorems of Prager 
and co-workers are essentially different from my theorems relating to the 
yield point of a plastic-rigid body, and that they have a greater practical 
value. My view is that the theorems are formally identical and give 
exactly the same results, 

The stress in a p-r body under specified surface tractions is obtained by 
indefinitely increasing the rigidity modulus G in a p-e solution. The 
influence of G is likely to be negligible whenever the plastic strain exceeds 
the elastic strain by a factor of order 3 (the straining being reasonably 
monotonic). However, these are by no means the only circumstances (as 
Dr. Lee maintains) where p-r analysis gives useful results. It has been 
proved, for instance, that in a wide class of plane strain problems 
the stresses are independent of G, even when plastic strains are of elastic 
order (Hill 1950 b, p. 244). This is true also in a bar under given torque 
or in a hollow sphere under given pressure (Hill 1950 b, pp. 87 and 99) 


and naa the consent of the author, Dr. Lee’s paper was shown by the editors to 
. Hill. 
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A precise meaning is attached to the limit or yield-point load of a p-r 
body: it corresponds to the moment when local distortion jirst begins 
under the given loading programme. By contrast, the load-distortion 
curve for a p-e body usually ascends smoothly from zero, and only when 
it passes through or attains a stationary value can a precise meaning be 
attached to the term p-e yield point. If changes of geometry and hard- 
ness are neglected in the p-e analysis, then, under reasonably monotonic 
loading, the load-distortion curve approaches the p-r yield point load 
asymptotically, though it does not usually reach it in any finite strain (e.g. 
Hill and Siebel, 1951, Phil. Mag., 42, 722). (The p-r work-rate is approached 
from below, according to the maximum work principle.) When geometry 
changes are not neglected, the p-e curve has a more or less pronounced 
knee to which the p-r yield point roughly corresponds. The p-r yield point 
for the undeformed body has, therefore, a clear practical significance. Dr. 
Lee perversely maintains, however, that the p-r analysis must be based 
on the deformed body after some unspecified strain, while allowing that the 
p-e analysis based on the wndeformed body is satisfactory. 

Consider, next, the maximum work principle and its complement (Hill 
1950 a and 1950b, p. 66). The associated continued inequality sets 
arbitrarily close bounds to the work-rate at the yield point of a p-r body. 
This is analytically identical with that afterwards given by Drucker, 
Greenberg and Prager (1951), and Prager and Hodge (1951, equations 
37.23, 39.12, 39.14). This later work, being restricted to proportional 
loading of non-hardening material, is on these two counts a special case of 
my 1950 theorems; the results of Greenberg and Prager (1951) for a 
structure of hinged bars are still more specialized. On the other hand, 
this later work was more general than mine in that it was formulated for a 
partly-plastic body ; examples of its use were also given. 

However, Prager stated, without reservation (e.g. Prager and Hodge 
1951), that the continued inequality provides arbitrarily close bounds to 
the p-e limit load, but only proved this 7f deformation occurs under 
constant load when geometry changes are neglected. As already remarked, 
such a stationary load is usually approached at infinite strain and then only 
if the straining is reasonably monotonic ; if boundary motion is allowed 
for, the curve rises steadily when the system is stable. Thus, Prager’s 
interpretation is approximate to a degree that is impossible to state 
beforehand. The only precise meaning that can invariably be given to the 
load concerned in these theorems is that it is the yield point of a p-r body. 
Accordingly I disagree entirely with Dr. Lee’s opinion that ‘the p-e 
approach provides definite information, whereas nothing can be deduced 
with certainty from the p-r theory’. Elsewhere (Drucker, Greenberg, 
‘Lee and Prager 1951) Dr. Lee seemingly contradicts himself: ‘the 
application of a complete p-r solution in the theorems gives the exact 
limit load for the p-e problem’. Furthermore, one of the illustrations 
cited was one where the load rises steadily (compression of a block between 


rough plates). 
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LI. Growth Spirals Originating from Screw Dislocations on Gold 
Crystals 


By 8S. AMELINCKX 
Geological Institute, University of Ghent* 


[Received January 11, 1952, revised MS. January 28, 1952] 


SUMMARY 


The present paper gives a more detailed account of the observation of 
growth spirals on (111) faces of gold crystals, obtained by precipitation. 
The observations prove that Frank’s theory also holds for metal crystals ; 
moreover they constitute the first visual proof for the presence in metal 
crystals of dislocations in general, and screw dislocations in particular. 
They further show that the Burgers vector of dislocations in metals is 
not necessarily limited to one unit of a few simple lattice displacements, as 
is generally assumed. The shape of the growth spirals is briefly discussed. 


§1. INTRODUCTION 


ACCORDING to the theory of Burton, Cabrera and Frank (1949, 1951) 
crystals can only grow, at supersaturations too low for two-dimensional 
nucleation, when dislocations are present. Where a screw dislocation 
ends in a crystal face an exposed edge is formed which, whilst growing, will 
subsist and wind itself up in the shape of a spiral (Frank 1949). Recently 
several experimental confirmations of this theory have been found on 
crystals of beryl (Griffin 1950), m-hexatriacontane (Dawson and Vand 
1951), carborundum (Amelinckx 1951 a, ¢; Verma 1951 a, b), and 
cadmium iodide (Forty 1951). In the case of mica the presence of dis- 
location lines with a screw component normal to the cleavage face has been 
proved by means of multiple beam interferometry (Amelinckx, in the press), 
Arguments strongly suggest that this is also the case for other layer lattice 
minerals such as kaolinite, nacrite, dickite (Dekeyser and Amelinckx 1951). 

Although the concept of a dislocation has been introduced to explain the 
mechanical properties of metals, no direct experimental evidence has been 
obtained so far for their existence in metals. The present paper gives a 
more detailed account of the observation of growth spirals on gold crystals, 
obtained by precipitation, and offers as a consequence an almost direct 
visual proof for the presence of screw dislocations in metal crystals. 
A preliminary account of these observations, in which the preparation of the 
gold crystals in question was described, has been published elsewhere 
(Amelinckx, Grosjean and Dekeyser, 1952). 


* Communicated by Professor W. Dekeyser. 
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§2. DEscriIpTION oF CRYSTALS 


The crystals were tabular ; they measured 5-40 microns across, and had 
hexagonal or trigonal shape. In fact they were octahedra developed in a 
tabular way following (111). Gold, having a face centred cubic lattice 
(a=4.07 A), can be considered as a close packing of spheres. The (111) 
planes are the close-packed and consequently the normal habit faces ; at 
low supersaturations such faces can only grow by means of spiral growth, 
As described in detail elsewhere (loc. cit.), the crystals were obtained by 
boiling a standard gold plating solution with concentrated hydrochloric 
acid. ‘When this mixture had cooled, crystals of sodium chloride and of 
gold were found on the bottom of the vessel, after a few days. It was 
therefore suggested by Prof. De Block (private communication), that the 
gold crystals had probably grown by epitaxy. It is indeed known (Van 
der Merwe 1949) that gold crystallizes with its (111) face on (110) of NaCl. 
The importance of this remark will appear afterwards (§ 6). 


§3. TECHNIQUE OF OBSERVATION 


The crystal faces were observed by means of an ordinary metallurgical 
microscope using bright field illumination. The photographs were taken 
slightly off focus, as this procedure increases the visibility of small surface 
steps. It has been shown that steps of only 7A high are still visible by 
this method, when conditions are good (Amelinckx 1951 b). As most of 
the crystals had their (111) face not exactly normal to the incident light 
beam, conditions were not favourable for observation. As a consequence 
only spirals with a relatively great step height could be photographed. 

Owing to the small dimensions of the crystals it was necessary to use the 
highest magnification. The consequent short working distance of the 
objective prohibited, unfortunately, all conventional interferometric 
measurements. The following procedure however enabled an estimate of 
the heights of steps. A small drop of a water-alcohol mixture was put on 
the crystal face. After a few seconds this drop evaporates partially, 
leaving a thin layer of liquid. Between the interfaces liquid-crystal and 
liquid-air two-beam interference fringes of the Fizeau type are formed,when 
using parallel, monochromatic light at normal incidence. Greatest 
dispersion is obtained just before complete evaporation. Unfortunately 
this procedure only allows visual interpretation as the fringes are constantly 


moving. 


$4, SHAPE AND ORIENTATION OF GROWTH SPIRALS 


On approximately 10°, of the gold crystals investigated, growth spirals 
with a well defined centre have been observed on the octahedral faces ; 
typical photographs are shown on plate XXIV.* 

In comparison with the growth spirals visible on carborundum their 
development is poor. In most cases only a few turns can be seen. The 


* For plate see end of issue. 
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spiral growth fronts are not smooth but ‘kinky’ (figs. 2, 3, pl. XXIV). The 
degree of ‘ polygonization ’ (Amelinckx 1951 ¢) is always very near to 
unity. 

The shape of the spiralsis dependent upon the rate of advance of a growth 
front in different crystallographic directions. It will therefore be triangular 
as a consequence of the triad symmetry axis perpendicular to (111) in 
cubic crystals. The loops of the spiral are parallel to the edges of the 
octahedral face (fig. 8 (a) and figs. 1,2, pl. XXIV). There are thus three main 
minima for the radial growth velocity corresponding to the directions 
perpendicular to the edges of the octahedral faces marked (a), (6) and (c) on 
fig. 9. These directions are the most close-packed. 


Fig. 8 


o 


Idealized shape and orientation of growth spirals. 
(a) triangular spiral, main minima only determine the shape. 
(b) intermediate shape : main and secondary minima are determining the 
spiral’s shape. 
(c) hexagonal spiral: the secondary minima dominate the main minima. 
(d) probable orientation of growth spiral on hypothetical hexagonally 
close-packed gold. ; 


There are six secondary minima (1, 2,..6). In certain circumstances, 
which are not clear for the moment, these secondary minima influence the 
shape of the spiral (figs. 5, 6, pl. XXIV). The symmetry remains trigonal ; 
the idealized shape becomes as represented on fig. 8(b). In the extreme 
case when the influence of the secondary minima is dominating, the 
shape becomes hexagonal (fig. 8 (c) and fig. 7, pl. XXIV); the hexagonal 
spiral being rotated 30° with respect to the hexagonal outline of the 
crystal. Polar diagrams for the radial growth velocity, as deduced from 
the shape of the spirals, are represented on fig. 10 for the two extreme 
cases (a) for fig. 1, pl. XXIV, and (bd) for fig. 7, pl. XXIV. 


In rare cases hexagonal spirals with their branches parallel to the face’s 
edges were observed (fig. 4, pl. XXIV). It is suggested that these could be 


. 
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hexagonally close-packed crystals, making the symmetry on (111) really 
hexagonal. It is indeed sufficient that the sequence ABAB .. . is formed 
once, by mistake, and that a suitable dislocation is present, to make sure 
that this arrangement will be repeated throughout the crystal or a part of 


Fig. 9 


Crystal structure of gold seen along [111]. The directions of main (a, b and c) 
and secondary (1, 2, 3,..6) minima in growth velocity are shown with 
respect to the lattice and to the outline of the octahedral face. 


Fig. 10 
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Polar diagrams for the radial growth velocity as a function of the crystallographic 
direction in the (111) plane, for the spirals of fig. 1 (a) and fig. 7 (b) of 
plate XXIV. 


it (Amelinckx 1951 c). This could however not be checked, the crystals 
being too small for single crystal x-ray diffraction. A powder diagram 
of the bulk material showed only the lines of the cubic lattice. 

The kind of mistake mentioned above is probably one of the origins of 
dislocations. 
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To estimate the heights of the steps the interferometric procedure 
described above was used. The deviation of the Fizeau-fringes caused by 
a known number of steps was estimated. In one case we found a step 
height of the order of 400 A--100 A (fig. 6, pl. XXIV). In other cases no de- 
viation of the fringe-system could be seen (figs. 1, 5, pl. XXIV) proving that 
the steps are less than 100 A. 

These measurements prove that the spirals cannot be due to dislocations 
of unit strength (d,,,—2, 2 A), but must be generated by groups of dis- 
locations all of one sign, having resultant Burgers vectors with a screw 
component, normal to (111), of a rather great number of repeat distances 
long. 

§5. TOPOLOGICAL CONSIDERATIONS 


In view of the calculation of the normal growth velocity (velocity normal 
to the crystal face) it is necessary to have information about the topological 
relations between dislocations, when the spirals generated by them are 
interacting. In the case of gold we only observed single or multiple 
growth promoting dislocations generating only one spiral dominating the 
growth of the face. Left- and right-handed spirals were seen in approxi- 
mately equal numbers. 

No spirals of opposite sign, interacting so as to form closed terraces, 
were seen, although special attention was paid to this and a very large 
number of crystals were investigated. A possible explanation is as follows. 
Dislocations of opposite sign attract each other and try to destroy one 
another, the attraction being roughly proportional to the product of the 
strengths of the dislocations and to the inverse of the distance. In the case. 
of carborundum dislocations of opposite sign, situated very near one to 
another, can be observed ; in metals however the mobility of dislocations 
is much greater, so that it is probable that dislocations of opposite sign and 
of suitable strength to produce visible spirals, cannot coexist within the 
small growing nucleus. 


§6. Distocation TUBES 


Very recently Frank (1951) has shown theoretically that when the 
Burgers vector of a screw dislocation exceeds a certain critical value, (10 A), 
the dislocation line is the axis of a hollow tube, in the interior of the crystal. 
The end of this tube has the shape of a crater. Frank considers the dots 
seen in the centres of growth spirals of carborundum as being caused by 
these craters, thus confirming his theory. He states however that these 
craters are not likely to occur in metals. It is now clearly visible on our 
photographs that the point of emergence of the dislocation line is marked by 
a dot. Apparently a discrepancy exists between theoretical expectation 
-and the interpretation of observations. The reason for this is simply 

that Frank did not consider the possibility of dislocations having a Burgers 
vector exceeding 10 A in metal crystals. Dislocations of the value 
mentioned above (400 A) are indeed unlikely to occur in metal crystals, 
unless they have grown in very special circumstances. Such circumstances 
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are for instance, according to Frank (private communication), growth on a 
foreign crystallization nucleus, in which case the gold crystal would 
inherit the screw dislocations present in the host crystal (this is in accord- 
ance with the growth process postulated in § 2). Once the large hollow 
screw dislocations are present, they are stable, as a large activation energy 
is necessary to detach unit dislocations from them. 

It seems to us that the large amount of misfit (25°) between the gold 
lattice and the host lattice is also responsible, at least partly, for the 
creation of dislocations. 
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LII. Baperiments on Two-dimensional Supersonic Flow in Corners and 
over Concave Surfaces 


By N. H. JoHANNESEN 
Fluid Motion Laboratory, University of Manchester * 


[Received February 4, 1952] 


SUMMARY 


Two simple cases of two-dimensional supersonic flow were investigated 
by schlieren photography, viz. the flow in two consecutive corners and 
the flow over a circular-are profile concave to the stream. The flows 
were photographed for smooth and rough surfaces of the models, Le. 
with laminar and turbulent boundary layers, and the observed flow 
patterns are compared with those predicted by inviscid theory. The 
flow with a turbulent boundary layer was in good agreement with inviscid 
theory outside the boundary layer, whereas the laminar boundary layer 
separated and a shock wave originated at the point of separation. In 
both cases the agreement was good at points far from the surface. 

A discussion is given of the flow at the point of intersection of three 
shock waves, and of a shock wave and a centred compression. For the 
latter case photographs show a weak shock wave in the direction of the 
downstream Mach line at the point of intersection. 


$1. INTRODUCTION 


WHEN a _ uniform two-dimensional supersonic stream meets two 
consecutive sharp corners, inviscid theory predicts two shock waves 
originating at the corners and combining into one shock wave at some 
distance from the surface. Further, since in general this simple three- 
shock configuration is not by itself a possible solution of the flow 
equations, it is necessary to introduce a weak expansion or shock wave at 
the point of intersection. 

When the uniform stream meets a circular-are profile which is concave 
to the flow and has its initial direction parallel to the free stream, inviscid 
theory predicts a shock-free flow (compression in a simple wave) near the 
surface, with a shock wave starting on the leading Mach line at some 
distance from the surface. 

These two cases of inviscid flow are considered in §3 below. The 
discussion includes a treatment of the flow near the point at which a 
vortex sheet originates as the result of the interaction of three shock-waves 
or of a shock wave and a centred compression. For the latter case the 
sign of the weak disturbance (compression or expansion) at the point of 


intersection is determined. The flow over the circular-are profile is 
ee ES ee ee eee 


* Communicated by W. A. Mair. 
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examined by combining the characteristic equations for rotational flow 
and the shock-wave equations. A numerical example is computed for 
a free-stream Mach number of 1.97. 

In §4 the experiments are described and the results compared with 
the inviscid flows. The flows were photographed for smooth and rough 
surfaces of the models, i.e. with laminar and turbulent boundary layers. 
The upstream effect of the shock wave or compression was found to be 
different in the two cases. In some cases the weak disturbance mentioned 
above could be detected. 


§2. Noration 
M Mach number. 
P ratio of local pressure to free-stream pressure. 
FR gas constant. 
a velocity of sound. 
p pressure. 
q velocity. 
ry radius of circular-are profile. 
Ss entropy. 
¢ function defined in eqn. (6). 
x, y Cartesian coordinates. 
y ratio of specific heats. 
6,, 55 deflections at first and second corners. 
# direction of stream. 
pw Mach angle. 


p density. 
es. 
ri sca 
2yR 


Suffix 1 denotes free-stream conditions. Other suffixes are defined in 
the figures. 
§3. THe Inviscip FLows 


We assume the free stream to be uniform and the flow over the obstacles 
to be adiabatic. The entropy is assumed to be constant along each 
streamline except through a shock wave. The effects of viscosity, and 
thus the boundary layers, are ignored. We deal with the two cases 
separately. 

(a) Flow in two consecutive corners. Let us consider the flow over 
a profile consisting of a small wedge placed on top of a larger wedge 
(fig. 1), and let the deflections be 6, and 6,. The two shock waves in the 
corners are geometrically fixed, ab the Pendens in region 3 depend on 
the free-stream conditions and on 8, and 8,. We shall assume that the 
flow in 3 is supersonic. 

The two shock waves combine into one, and a vortex sheet originates 
at the point of intersection. Across the vortex sheet the pressure and 
flow direction must be continuous. These conditions can, in general, 
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not be satisfied by such a three-shock configuration, and it is necessary 
to introduce an additional weak shock wave or expansion at the point 
of intersection. This weak disturbance must lie between the shock 
wave (2-3) and the vortex sheet, separating the regions 3 and 4. Its 
direction is approximately that of a Mach line. 


> = 
<SS ONS = 
Le 


Flow over double wedge. 


Let p denote the pressure, P the ratio of pressure to free-stream 
pressure p,,6 the velocity direction measured from the free-stream 
direction, and M the Mach number. For the single oblique shock wave 
(1-5) P; can be plotted against 6, for a fixed WM, (fig. 2). The point 
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Relations between pressure ratios and streamline deflections. 


A with coordinates (03, Ps) will, in general, not lie on the curve P;=P;(@;). 
The points (6,,P,) for various strengths of the weak disturbance 
(3-4) will then lie on the curve BAC, where the branch 4B corresponds 
to an expansion and AC to a compression. The position of A in relation 
to the curve P;—P,(8;) determines whether the disturbance is an 
expansion or a compression. The small deflection at the disturbance 
is 8p —04, where D is the point of intersection of the two curves. 
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For the values of M,, 5, and 5, of the experiments described below the 
disturbance was found to be a compression, but its strength was so small 
that the shock-wave could not be seen on the photographs. Thusfor M eh 
6,=5° and 8,=10°, 0)—0,= —0-02° and (p,—;)/p,=0-001. This is 
only about a tenth of the minimum strength that can be detected by the 
schlieren system. 

The theory of the three-shock configuration was not investigated in 
detail, but the simpler case where the two corner shock-waves are replaced 
by a continuous centred compression is considered in the following 
section. 

(6) Flow over a concave surface. Assume a uniform stream of Mach 
number M, meeting a circular-are profile which is concave to the flow 
and has its initial direction parallel to the free stream (fig. 3). Near the 


iors 


Flow over concave surface. 


surface the compression takes place through a simple wave, but the Mach 
lines from the surface will converge and intersect, and a shock wave will 
be formed at some distance from the surface. 

The flow pattern can be computed by combining the characteristic 
equations for rotational flow and the shock-wave equations. In two- 
dimensional rotational flow there are three families of characteristics : 
the plus and minus Mach lines* and the streamlines. 

With Cartesian coordinates x,y the characteristic equations for the 
Mach lines can be written (see Meyer (1948 a), pp. 205-206) 


dy ee n8) 
dg tan (OTH); Sco Ce Oh) aCe (1) 


sin ps COS 
pa? 


do dp=0, : oPbro Fol eames. (2) 


* An observer looking downstream from a given point has the plus Mach line 
through the point on his right and the minus Mach line on his left. 
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where the upper and lower signs refer to a plus and a minus Mach line 
respectively. § is the angle between the z-axis and the velocity direction, 
u. the Mach angle, p the density, p the pressure and a the local velocity 
of sound. 

For the streamlines we have 


oy ow ds=0, CE eS eee 
du 

2d 
qdq+ —* =0, Jt, 2 ae 


where s is the entropy and q the velocity magnitude. 
Equation (2) can be written in the form 


d(t+0)-+do sin 2u=0 (5) 


with dt=(dq/q) cot . and do=(ds)/(2yR), where y is the ratio of the specific 
heats and # the gas constant. 
It is found that 


t=pu—n/2+dA-*! cot! (A~! tan p) PY one eck LOD 


with A=[(y—1)/(y+1)]}'?. The function (1 000—#), with ¢ in degrees, 
has been tabulated (see Herbert and Older (1946)). 

The shock-wave equations are most conveniently written so that the 
conditions downstream of the shock are expressed in terms of the shock- 
wave angle and the conditions upstream of the shock wave (e.g. see 
Ferri (1949), p. 45). 

When the flow is isentropic, the streamlines are no longer characteristics. 
The equations of the Mach lines become 

d(t--0)=0, - =tan (0-Fu). oir ooryyr) 4 i 4: 

Near the curved surface the flow is isentropic. In fig. 3 let OF be the 
leading Mach line, FGH a minus Mach line, and KF and HG@ the plus 
Mach lines through F and G. Using (7) we find 


ty +0,=tp+Op=tgtOg, 


aoa (| 
tg—Og=tp—Op, | oy 
or O,—6q and tp=tg, i.e. wp=ug, Which shows that near the wall the 
minus Mach lines are straight. These Mach lines form an envelope as 
shown in fig. 4. The envelope starts at a point P on the leading Mach 
line. The position of P is found as the limiting position of # when the 
Mach line FE approaches the leading Mach line OZ. It is found that 

eS es 

Yp= y+ e “Me A Bees, he Oe ena te aa te (9) 


where r is the radius of the circular arc. This result agrees with the 
condition for shock-free flow at the inlet of a two-dimensional contraction, 
derived by Meyer (1948 b, p. 455). 
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The envelope is situated upstream of the leading Mach line, which 
may be considered the second branch of the envelope (see Courant and 
Friedrichs (1948), p. 115). In order to obtain a simple approximate 
method of predicting the shock position, the ratio a/b, fig. 4, defining the 
position of the shock wave (shown as a broken line) in relation to the two 
branches of the envelope, was determined under the following assumptions: 
(i) the Mach lines are straight, (ii) the slope of the shock wave is the mean 
of the slopes of the Mach lines upstream and downstream of the shock 
wave, and (iii) the shock wave has the shape y, «<,% (see the x,, y, system 
in fig. 4). The ratio a/b was found to be 0-75 for «=2, 0-71 for «=3/2, 
and 0-78 for «=. The value of 0-75 would therefore be expected to be 
nearly correct under widely different conditions. 


Fig. 4 


Envelope of Mach lines. 


The flow behind the shock wave is rotational, but the entropy is 
constant along each streamline, and the flow can therefore be computed 
by using eqns. (3)-(6). In the general case of a supersonic rotational 
flow, it is necessary to cover the flow field with a network of Mach Jines 
and streamlines and compute the properties of the flow and the geometry 
of the network simultaneously (see Ferri (1949), p. 86). However, in our 
case it is possible to use an easier procedure by which the properties at 
all the points in the network can be determined without knowing the 
positions of these: points, i.e. the geometry. The geometry can then be 
determined afterwards by simple trigonometry. 

In fig. 5 let PQR be the shock wave, REMSN a minus Mach line, PS the 
plus Mach line through the point P (at which the shock wave starts), 
PM the streamline through P, and QL a streamline between P and R. 
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We then have from (5) with suffix R denoting conditions downstream 

of the shock wave at R 
and from (8) 

05=6y ; ty=ty >; Msn 5 ts +0 s=t,. . . 5 ° (11) 
Let the free-stream values s, and o, be zero. As a first approximation 
we then have 

[do sin 2u]§=op sin (up+s) 

and therefore 

tp—Op—2tg tt, +op SM (URt+Hs)=9. sey ace a2) 
For known conditions at R, ws and ty can be determined by trial and 
error, using (12) and (6). 6, is found from (11). 


Fig. 5 


Computation of flow pattern. 


For any point J on the minus Mach line between Rand S we have 
(tp—Op)—(t, —0,) +[do sin 2u]}'=—0, 


(t,—01,) —(ts—O5) +[do sin 2n]s=0 
or 


and, approximately 
[do sin 2\7=(ep—2) sin (up+ez), 
[do sin 2u]§=ag sin (uz+ps)- 


Eqn. (13) gives 4, by trial and error, and the remaining conditions at 
ZL, can then be found. 


di : : 
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When sin (up+4z) and sin (uz+ys) are nearly equal to unity, 4, cannot 
be determined accurately from eqn. (13), but must be found by 
interpolation from the values at neighbouring points. 

Fig. 6 shows the flow that was computed for M,=1-97, y=7 /5 and a total 
stream-deflection of 22-5° at the curved surface. (At this value of 
M,, the velocity would be sonic after a deflection of 25-5°.) The ratio 
a/b as defined above was found to be 0-75 to 0-80, in good agreement with 
the simple estimate. 

Because of the big steps used in the computations, fig. 6 could only be 
expected to show the general character of the flow pattern predicted by 
inviscid theory. However, the construction of the flow pattern was 


Fig. 6 


\ 


Flow over circular-arc profile for M,=1-97. 


continued graphically, assuming that there is an expansion region with 
centre at 7' (fig. 6), and when the position of the decaying shock wave 
was compared with the corresponding experiments, the agreement was 
found to be good (see below). 

In order to investigate the possibility of a flow with a detectable weak 
shock wave at the point where a vortex sheet starts, we consider the case 
where the curved surface is so shaped that the continuous compression 
is centred and the envelope reduced to a point (fig. 7). Such a flow is 
physically possible, and can be obtained by giving the fixed boundary 
the shape of a streamline in a suitable Prandtl—Meyer expansion. 
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With the notation of fig. 7, we find the curves P5=(p5/P1)=P5(95) and 
P,=P;(63) similar to the curve P,=P;(0;) in fig. 2. If we expand the 
pressure ratios P, and P,; in powers of the deflection 6, we get : 


Py=1+C,0+C2+ 0,0 +0(6), | (14) 
Ps=140,0+0,0?+ (C;—D)®+0(), 
where 
D= 1 yM (M21)? [(5—3y) Mt (12—4y) MP +8) - 5) 


(see Ferri (1949), p. 76, and Kahane and Lees (1948)*), D=0 for 


Fig. 7 


= 


Flow with centred compression. 


M,=1-245 and M,—2-540 and is negative between these two values of 
M,. Therefore, for small values of 6, the weak disturbance is a compression 
for 1:245<M,<2-540 and an expansion for M,<1-245 or M,>2-540. 

The exact values of M, for which (P;—P;)=0 depend on 6, and for 
high values of 6 they are slightly different from those given above. 

For the case of M,=1-97 and a deflection of 22-5° through the 
continuous compression (corresponding to the experiments described 
below) the weak shock wave was found to have a_ strength 
(P4—P3)/P3—90-08, and we should therefore expect to be able to see it on 
the schlieren photographs. 


/ 


* It should be noted that eqn. (15) is given incorrectly in Busemann’s origina | 
paper (Lufo. 12, 210-220 (1935)) and in a number of other papers and books. 
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$4. EXPERIMENTS 


The experiments were made in a small intermittent supersonic wind 
tunnel which has been described by Bardsley and Mair (1951). The 
working section of this tunnel is 12-510 em, the Mach number 1-96 and 
the Reynolds number 1-3 10° per cm. The humidity of the air was in 
all cases less than one part in 2 000 by weight and usually much lower. 
The models were made of high-carbon steel, ground and lapped to the 
correct shape. They extended over the whole width of the tunnel, and 
were supported by two arms on the downstream side. The position of 
these arms could be altered to give different angles of incidence of the 
models. The schlieren system was of the conventional two-mirror type, 
and a flash tube was used as light source because it was conveniently 
available. The knife edges were in all cases parallel to the free stream. 
The flows were also observed on a screen using a continuous light source, 
and were found to be steady. 

(a) The double wedge had the dimensions shown in fig. 8. Experiments 
with smooth surfaces were made for 6, approximately 0, 24, 5, 74 and 
10 degrees. Figures 9 and 10 (Plate X X V*) show the flows for 6,=—0-1° 
and 6,=7:2°. 5, was 9-4° in both cases. For 6,—10° the tunnel was choked: 


Fig. 8 


Section of double wedge (with dimensions in mm). 


Tn figs. 9 and 10 (Plate X XV) it can be seen that, in addition to the shock 
waves at the leading edge and in the corner, there is a shock wave from the 
point where the laminar boundary layer separates. The flow downstream 
of the separation shock-wave is apparently uniform and the edge of the 
boundary layer is therefore nearly straight. The dark line in the white 
boundary layer downstream of the separation point is due to ‘ over- 
loading’ of the schlieren system (see Mair (1952)) and could be made 
to disappear by increasing the width of the second slit. It should be noted, 
however, that it may be advantageous to adjust the schlieren system so 
that the dark line occurs, as it sometimes makes it easier to interpret the 
photographs. The overloading indicates a strong density gradient in 
a layer separated from the surface. In fig. 10 the vortex sheet starting 
at the point of intersection of the two strong shock waves can be seen. 


*For plates see end of issue. 
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By measuring the photographs, it was found that the slope of the 
first shock wave agreed with inviscid theory within 3°. The second 
shock wave, where the boundary layer separated, gave a deflection of 
approximately 5°. The position of the point of separation was 6 to 
7 mm upstream of the corner. This is approximately 60 times the 
displacement thickness of a boundary layer on a flat plate at the Mach 
number and Reynolds number of the experiments and at the same 
distance from the leading edge as the corner. The angle of the main 
shock (1-5) agreed within 0-2-0-5° with the value calculated for a total 
deflection of (6,-+6,). 

The surface was then made rough by spraying it with ‘dead black’ paint, 
and the photographs figs. 11 and 12 (Plate X XV) were obtained for 6, =—0-0° 
and §,—7-5° respectively, with 5, equal to 9-4° in both cases. In fig. 12 the 
upstream effect of the corner shock-wave is hardly detectable. In 
fig. 11 there is still some upstream effect, but much less than in the case 
of the smooth surface. It is believed that the difference in magnitude of 
the upstream effects shown in figs. 11 and 12 is not due to the difference 
in Mach number upstream of the corner shock-wave, but is caused by 


Fig. 13 


Section of circular-are plate (with dimensions in mm). 


the greater roughness of the surface in fig. 12. The method used for 
roughening the surface had the disadvantage that the paint was blown 
off in the course of a few seconds, and it was therefore not possible to 
get exactly the same roughness in the two cases. In both figures the 
wavelets from the rough surface are clearly visible. The irregular wavy 
line below the model in fig. 12 is a water smear running downstream 
along the glass window, and it therefore represents a streamline in the 
boundary layer on the window. Note that the trace bends before reaching 
the image of the shock wave from the leading edge, because the shock 
wave bends upstream when it runs into the boundary layer. The water 
smear was caused by a drop of dirty water issuing from a cavity in the 
upstream joint of the glass window. The water drop remained there 
after the window had been washed and dried, but it could only have 
a momentary effect on the dryness of the air in the main stream. 

The circular-are profile had the dimensions shown in fig. 13. A flat 
part upstream of the curved part prevented the shock wave due to the 
slight bluntness of the leading edge (see Bardsley (1951)) from obscuring 
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the flow over the curved profile. The total deflection of a streamline 
along the surface of the profile was 24-5°. With a smooth surface the 
flow was photographed with the straight part set at incidences 
0-2°, 2-7° and 5-1°. These incidences correspond to Mach numbers just 
upstream of the leading Mach line of the curved surface (i.e. in the 
effective free stream) of 1-97, 2-05 and 2-15 

Figures 14 and 15 (Plate XXVI) show che flow for M,=1-97 and fig. 16 
(Plate XXVI) the flow for M,=2-15. As on the double wedge the 
compression is propagated upstream, the boundary layer separates, and a 
shock wave occurs at the point of separation. Again the flow downstream 
of the separation shock-wave is apparently uniform, and the edge of the 
boundary layer is therefore nearly straight. 

By measuring the photographs, it was found that the shock-wave 
angles agreed well with the values calculated from the measured slopes 
of the boundary layer. The total deflection of the streamlines (outside 
the boundary layer) was found to be 22-5°, ie. 2° less than that of a 
streamline at the surface. This is because the expansion at the downstream 
edge of the model produces a negative pressure gradient in the boundary 
layer on the upstream side, so reducing the rate of growth of the boundary 
layer. For the decaying shock wave it was found that the distance 
between the computed wave and the one found experimentally* was 
only about 1 mm on the scale of the experiments. This corresponds to 
about 1/50 of the radius of the circular are. 

The deflection at the point of separation was approximately 4°, and 
there was a further deflection of about 4° at the point of reattachment. 
The point of separation was about 2 mm upstream of the beginning of 
the circular arc. This is about 20 times the displacement thickness of an 
undisturbed boundary layer at the beginning of the curved part. 

When the surface was made rough the photograph fig. 17 (Plate XX VI) 
was obtained for /,—1-96 (zero incidence). The upstream effect is no 
longer detectable and the flow everywhere outside the boundary layer is 
in good agreement with inviscid theory (fig. 6). 

In all the photographs figs. 14-17 the flow corresponds roughly to the 
case of a centred compression meeting a single strong shock wave, and the 
weak shock wave in the direction of the downstream Mach line can be 
seen in figs. 14 and 15. It could also be seen on the negatives from which 
figs. 16 and 17 were reproduced. 


§5. CONCLUSIONS 


The supersonic flow over a double wedge (fig. 8) and a circular-are 
profile (fig. 13) was photographed and the actual flow compared with 
that predicted by inviscid theory. 

It was found that when the surface of the profile was rough and the 
boundary layer therefore turbulent (e.g. see figs. 12 and 17), the flow 
outside the boundary layer was in good agreement with the inviscid flow. 


a 


* The downstream side of the image represents the true position of the shock 
wave outside the boundary layer (see Bardsley and Mair (1951)). 
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When the surface was smooth and the boundary layer laminar the increase 
of pressure on the profile was propagated upstream through the boundary 
layer and caused separation (figs. 10 and 14). A shock wave originated at 
the point of separation and caused considerable deviation from the inviscid 
flow. The flow at points far from the surface was in all cases in good 
agreement with the inviscid flow, whatever the state of the boundary 


layer. 
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SUMMARY 


An attempt is made in this paper to improve the theory of two electrons 
of opposite spin in a field of spherical symmetry such as exists in some 
atoms. The method follows closely the theory of paired electrons in 
a chemical bond, given recently by the authors. In order to take fully 
into account the spatial correlation of the electrons due to their 
electrostatic repulsion, the wave function is expressed as a series. Each 
member of the series contains a function of the radial distances of the 
electrons multiplied by a Legendre function of the angle which the 
electrons subtend at the nucleus. It is then found possible to obtain 
a series of coupled equations for the radial functions in which the angular 
coordinate does not appear. The exact theory is used to examine some 
approximate expressions which have been used for the wave functions 
of the helium atom. 

Spatial correlation may occur in two ways. It may be of the ‘ in-out’ 
type, one electron tending to be outside the other, or it may be of the 
angular type, the two electrons tending to be on opposite sides of the 
nucleus. The spatial correlation of the electrons in helium in its ground 
state is examined from this point of view. 


§1. INTRODUCTION 


One of the greatest difficulties in the theory of atomic and molecular 
structure is to take adequately into account. the electrostatic repulsion 
of the electrons. Because of this interaction it is not possible to describe 
satisfactorily the properties of many electron systems in ordinary three- 
dimensional space. Even a two-electron system such as a helium atom 
in its normal state has to be described in six dimensions (or eight 
dimensions, if the spin properties are included). This is because there 
is a correlation between the behaviour of the electrons and in order to 
describe the system fully, it is necessary to know the probability of finding 
simultaneously one electron in one element of space and a second electron 
in another prescribed element of space. It has not as yet been possible 
to devise an experiment which would disclose this inner property of an 
atom. All that diffraction experiments, for example, can demonstrate 
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is the average total distribution in three dimensions. This total electron 
densitv of a helium atom is known to be spherically symmetrical. ‘The 
hands of a clock provide a simple analogy. Each hand has on the average 
the same probability of being in a prescribed element of angle, and if a 
photograph were taken with a long exposure the image of the two hands 
would have circular symmetry. But there is another quantity which 
could be deduced about the hands, and that is the probability of finding 
them at any prescribed angle to each other. So in an atom or a molecule 
there is a probability of finding the electrons in given configurations 
relative to each other. : 

The correlation of electrons in molecules is of more direct interest than 
in atoms because it is a major factor in determining molecular shape. It 
is for this reason a subject of great theoretical interest. It is accordingly 
instructive to examine the spatial correlation of electrons in atoms if 
only as a guide to the more complex systems which exist in molecules. 

While the spatial correlation of electrons in atoms is not itself directly 
observable, it none the less affects other observable properties. The 
energy content of atoms depends directly on the average distance apart 
of the electrons. Attempts to calculate atomic energies involve implicitly, 
if not explicitly, methods of estimating, or of guessing, the distribution 
of electrons relative to each other. Many attempts have been made to 
calculate the energy of two electron systems, like helium, by variational 
methods, but the primary emphasis has been on such observables as the 
ionization potential and not on the distribution of the electrons relative 
to each other. 

In this paper the effect of the electrostatic repulsion on the wave 
function of a. pair of electrons in an atom, when this has a spherically 
symmetrical field, has been examined by a method similar to that used 
elsewhere for a pair of electrons in a chemical bond (Lennard-Jones and 
Pople 1951). The wave function is expanded in a series similar to that 
used for the electrostatic repulsion when expressed in terms of the radial 
coordinates of the two electrons and the angle they subtend at the 
nucleus. The dependence on this angle then appears in terms of spherical 
harmonics, and it is necessary only to find how the wave function depends 
on the radial coordinates. It is shown how in principle an exact solution 
of this problem can be obtained. This theory is then used to examine 
the limitations of various approximate methods which have been used, 
and in particular to show the significance of the calculations of Taylor 
and Parr of the energy of the helium atom (Taylor and Parr 1952). 

For paired electrons in the lowest energy state of a system like helium 
it is the electrostatic repulsion which is responsible for spatial correlation. 
This is always operative whether electrons have the same spin or opposite 
spin. When electrons have the same spin, there is another factor which 
may be more important, and that is the exclusion principle. Though not 
of direct concern in this paper, the great influence of the exclusion 
principle on the spatial correlation of electrons must for the sake of 
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completeness be pointed out. This principle is usually expressed 
mathematically in the form that the total wave function of a many- 
electron system must be antisymmetric for interchanges of electrons. 
The significance of the principle is most easily appreciated by considering 
the fictitious system of many electrons moving without mutual repulsion 
in the same field, whether atomic or molecular. The wave function of 
this system can then be expressed accurately as a product, or sum of 
products, of one electron functions. For the lowest state the wave 
function has the form of a determinant. The exclusion principle then 
demands that not more than one electron of the same spin be assigned to 
the same one electron function. The wave function vanishes when any 
two electrons have the same spin and the same spatial coordinates. This 
means that electrons of the same spin avoid each other and tend to occupy 
different regions of space. The significance of this has been examined 
for a number of atoms and molecules in another paper (Lennard Jones 
1952). 


§2. THE SEPARATION OF THE WAVE. EQuATION 


The Hamiltonian of a system of two electrons moving in a spherically 
symmetric field V is, in atomic units, 


#=H,4+4H8.4+1/r4., . . e ° ‘ a ‘ (2.1) 
where 


H,=— 2+ V(r), | (2.2) 


A= pear V(r). 


In these formulae 7, and r, are the distances of the two electrons from 
the centre of the symmetrical field and r,, is the interelectronic distance, 
V may represent the field due to a central nucleus, as in the helium atom 
or the Lit ion, or it may represent the combined effect of a nucleus and 
the screening by a closed shell of electrons, as in atoms such as beryllium 
and magnesium. 

The interelectronic potential (1/7,.) depends only on 7,, 7, and 6,5, the 
angle subtended at the centre by the two electrons. It may be expanded 
in terms of Legendre functions in the form 


Wine) hh 4) (e08,0;5); Pr alia / 4) (205) 
1=0 
where r_ and r. are the smaller and greater of 7, and r, respectively, 
and where P, are Legendre polynomials defined so that P)(1)=1. For 
brevity this will be written 


1/ty3= U (74, 72)P 7 (C08 O45), ‘aurea are) 


l 


ius 


U , being independent of 65. 
ZR2 
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The solutions of the wave equation, when the Hamiltonian has the 
form (2.1), will be of various symmetry types. For those of symmetry 
18 it is convenient to expand y in a series similar to (2.4). Thus we write 


HL, 2)= obs %g)P (C08 O35); 2. hee) 
=0 


ys, being a function of 7, and 7, only. 

The expansions (2.4) and (2.5) are independent of the choice of a fixed 
reference frame for specifying angular position, but they can easily be 
transformed to expressions involving surface harmonics referred to such 
a frame. If the spherical polar coordinates of the two electrons are 
(7,,9,, 6,) and (79, 09, dy) respectively, then* 

P (cos 6,2)=P, (cos 0,)P; (cos A) 


1 (l—m)! 
DF 
ee eae 


Using (2.5) the last term in the Schrédinger equation becomes 


P|" (cos 6,)P;'" (cos 65) cos m(d,—¢}). (2.6) 


(1/72) p= oe U PP (COS 81 5)Py, (COS 49). dist): Rae 


A product of Legendre polynomials can be expanded as a sum by the 
formula (Adams 1878) ; 


Min(d, m) , 


P(w)Pn (u)= 2 Ga reer) 9 - C - (2.8) 


r= 


where the ¢,,,‘") are numerical coefficients given by 


(2.9) 
and 
1305 ee (2e—1) 


= (2.10) 


Using this formula, the complete Schrédinger equation can now be 
separated into a set of coupled equations for the functions ¥,, since the 
coefficients of the various Legendre functions must separately vanish. 
Thus we get 


Lied 220) 1 +d 0 i+-1)f 1 1 
Sy Pa [ee a2 eee ak 
{ E ry (r a) x re? Ors ¢ a) 7: 2 fe 2p 2 


ae V(ry)+ Vir)| b+ 2’ oy ee tras Lor, . 4 (2.11) 


Role 


where the ¢,,,,") are defined to be zero if r is negative or non-integral. 
ee a ae ee eee 

* For states of the system which do not possess spherical symmetry it will be 
necessary to use (2.4) together with (2.6) and to expand the wave function in a 
form similar to (2.6), each term in the expansion being multiplied by a function 
ofr, and rz. It is hoped to consider this problem in a later paper, 
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We note, in particular, that the equation for y, takes the form 


R laelaso a0 Lid oer. 
=i bere 7) beaeelle lt +Vr)+ Vea) | 


eee Ly 
SiS en In bn bo: Ae vag Tees he PIS she (221-2) 
Thus the purely radial function % depends on the coefficients of the 
higher angular terms in the full expression for *(1, 2). 

Up to this point the treatment is exact except for any approximation 
implicit in the use of a smoothed potential function V. The eqns. (2.11) 
for 4, may not be soluble explicitly, but as in the case of a pair of electrons 
in a chemical bond (Lennard-Jones and Pople 1951), the accurate theory 
forms a useful background against which to examine the limitations of 
approximate wave functions. 


§3. APPROXIMATE SOLUTIONS IN TERMS OF ONE-ELECTRON FUNCTIONS 

In order to examine the physical significance of approximate wave 
functions, we shall consider first those which are independent of the 
angle @,,. This implies neglecting all but the first term %, in the 
expansion (2.5). The simplest approximation to %, with the necessary 
symmetry properties is a product of one-electron functions, 


DTA aa fo) aes Ge Rh SUDA ths ay asled Ov) 
If fy is chosen so as to minimize the energy integral [%hxy dr / [ype dz, 
then this is the simple Hartree-Fock function. According to this function 
each electron moves in the averaged field of the other, but- otherwise 
there is no correlation between them. No allowance is made for their 
tendency to avoid one another because of their electrostatic repulsion. 
The function 4 cannot, of course, allow for angular correlation, that 
is the tendency for the angle 6,, to be as large as possible. More accurate 
forms than (3.1) do, however, introduce radial correlation, the tendency 
for the two electrons to be at different radii. As a better approximation 
we may accordingly replace (3.1) by 
b=foal1)fo,a("2) —So,0("1)So,0("2)- eae 64) 
The negative sign is used in (3.2) as it is found later that this avoids the 
use of imaginary functions. As we are only interested in the function %p, 
we shall drop the zero suffix for the rest of this section. We then have 


fyebdr=S,,2+5,,2—2892, « . - - + (3.3) 
and fbb dr=2(SugH aa +S pull np —2S ant an) 
+{(aa|U)|aa)+(bb|U|6b)—2(aa|U,|5b)}, - . (3-4) 
where we have written 
SaaS ale | 
Hei A aan, a 235) 
(aa |Uy|66)=Sf fa fal2)Uofo(1)fo(2) da, dx, | 


and other similar expressions. 
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To obtain. the optimum forms for f, and f, we have to minimize the 
integral [Wi dr, subject to the condition 
[dab drs ae Se 
This leads to the following pair of coupled integro-differential equations 


eA Be —Suofos i { ( Up)aada rf ( Uoavf ot 
a, ee ane (ES,,—H av) So 


' meer eee: 
ALS 5fo—Sarslat + {(U5 Veple Oo)astat 
Je at fot+(E8ar—Aar)lar J 
where we have written 
Jed =Jf,(2 )Uofal2) dxs, . . . . . (3.8) 


and other similar expressions. The solution of (3.7) must be scaled so 
as to satisfy the normalization condition 


Sua? tS oy°—28an7=1. . . . . . > (3.9) 
We note that if we put /,—0 in the first of eqns. (3.7), the ordinary 
Hartree-Fock equation is obtained. There are two distinct 


approximations, therefore, in this function, one being the expression 
of 4 as a simple product in the form (3.1) and the other the neglect of all 
the terms in (2.12) other than yp. 

An alternative form of the function (3.2) can be obtained by 
transforming to new one-electron orbitals given by 


X =(fatto)/ V2; 
x’ =(fa—fo)/-V'2; 
so that the wave function becomes 


b= (1) x42) x’ (D2) 2 ae ase 


To obtain further 2 aed it is necessary to include approxi- 
mations to the functions , for />1. Such approximations may be 
simple products as in (3.1), giving a total wave function 


(3.10) 


b= 2 flridfilrs)P (008 612), ~\ tcp reat Spel Saket 
or sums of two simple products as in (3.2) 
p= etal al )F1,a("2) Si, ol? )S1,0(%2)} P (C08 09). - (3.13) 


Kquations for the optimum forms of any of these functions can be 
obtained by minimizing the integral [iy dr subject to the normalization 
condition, using the corresponding form for ys. It is to be noted that the 


wave function (3.13) can be written as a generalization of (3.11) in the 
form 


p= PLxAl) (2)+x,'(1 (2)}P,( (cos 6,5). Crest (sel) 
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§4. IntustrativeE NumMericaL Resutts roR THE HeLium AToM 


In order to illustrate the way in which the introduction of the higher 
components of the wave function lowers the total calculated energy by 
allowing for the tendency of the electrons to avoid one another, it is 
instructive to examine some results obtained recently for the helium 
atom by Russell Taylor and Parr (1952).* Their calculations are based 
on hydrogen-like atomic functions instead of on the one-electron orbitals 
discussed above. By the use of such functions the calculations can be 
done simply, but it should be emphasized that there is no reason to 
suppose that they are good approximations to the optimum forms. 
The functions considered are the following :— 


Function 1. An approximation to %(r,,7.) of the form f,(1)fo(2) 
using a hydrogen-like (1s) function (with the best value of Z) for So: The 
optimum value of Z is 1.6875. 


Function 2. i: approximation to 4(7;,72) of the form 
{xo(1)x0'(2) + x0/(1)xo(2)} using hydrogen-like (1s) functions, with different 
screening eee , and Z’) for x, and x,’._ The values of Z and Z’ are 
1.19 and 2.18365. 

Function 3. To function 2 is added an expression f,(1)/,(2)P (cos 0,9) 
using a function proportional to the radial part of a ate like 
2p-function for f,. Expanding P,(cos @,,) by the addition theorem for 
Legendre polynomials, it is seen that function 3 corresponds to the 
superposition of the configurations (ls, 1s’) and (2p)?. The screening 
constant in the radial function /,, is taken to be 4.9504. 

Function 4. Similar approximate forms are added to represent the 
effect of (3d) and (4f) functions. The full expression used can be written 
‘in the form 

$={xo0(1) x0’ ae (2)} +A) ee 91) 
+-fe(1)f2(2)P2(cos paca is g(COSA,5). . . (4.1) 
The functions /, and f, are oe to be Ee to the radial parts 
of hydrogen-like (3d) and (4f) functions. The screening constants in 
the radial functions f, and f, are taken to be 11.305 and 19.04. 

The energies of the helium atom calculated from these functions are 

shown in table 1, which is taken from the paper of Taylor and Parr (1952). 


Table 1. Calculated Energies for the Normal Helium Atom 


— Ei calc Beale —EHexpt 
(units of e?/2a,) (units of e?/2a,) 
Function 1. 569531 0-1113 
Function 2. 5.751316 0-0552 
Function 3. 5:790369 0-0162 
Function 4. — §:794860 0-0117 


* We are indebted to the authors for allowing us to refer to their results 
before publication. 
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As Taylor and Parr pointed out, these figures illustrate the relative 
improvements obtained by allowing for the various types of correlation. 
The analysis of this paper shows the relation of these improvements to 
the properties of the exact wave function (2.5). 


$5. Tue Sparta, CORRELATION OF ELECTRONS IN THE HELIUM ATOM 


Some further insight into the correlation between the electrons in the 
helium atom can be obtained by a more detailed examination of the 
various parts of the wave function. It should be emphasized, however, 
that any conclusions based on the approximate functions of Taylor and 
Parr may be inaccurate because of the particular types of radial functions 
they have used. 

In terms of the accurate function (2.5), the probability of the radial 
distance of the first electron being between r, and 7,-+-dr,, of the radial 
distance of the second electron being between r, and 7,+dr, and of the 
angle subtended at the nucleus being between 6,, and 6,,+-d@,, is 


P(rs, 5; 035) dry drs db, 5, 
where 


P(11, To, 942)= I XP {P( cos 1) }? 
1=0 


+2°2 Pmibnl m(COS 8; 2)P,,(COs a.) | 87°r,?ro? sin 64. (5.1) 
m>n 

To investigate radial correlation irrespective of angular correlation 

in the helium atom, we must integrate (5.1) over all values of 0,, between 

0 and 7. This leads to a radial probability distribution function 
P(r,, ’2) dr, dr, where 


P(r o's)=167 ae 


a ee (5.2) 


From the form of (5.2) it is clear that the radial correlation will be mainly 
described by the first term of the wave function yp, partly because this 
is the largest term in the expansion (2.5) and partly because of the 
numerical factors (2/+-1)~?. 

If we neglect all but the first term in (5.2) and then use the particular 
approximation to % given by function 2 of the last paragraph where 


Zr } 


oc eo ) 
0) ahd ee eee 
Xo ce“, J 
then P(r,, 7) is proportional to the function 
F(s, t)=4(s?—#?)? e- "+7 cosh {(Z'—Z)t}}?, . . (5.4) 
where . : 
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It is of interest to find out under what circumstances the values of 
r, and r, tend to differ. From (5.4) we obtain 


2k 
oe |, =8s"{ (Z' —Z)?s?—2} eT *™2t2 (5.6) 


and this gives the sign of the curvature of the surface F along the section 
r,=Pr,. If s is large, this will be positive. In other words, if the average 
of the two radii is large, the two electrons will tend to separate into 
different spherical shells and radial correlation of the ‘in-out’ type 
will occur. But if s is small, (5.6) will be negative and the electrons will 
most probably be found in the same spherical shell. The critical value 
separating the two cases is 


ae, 2) Zier |e eet ee tie dues ee 8 (D7) 
Using the values of Z’ and Z given by Taylor and Parr, this corresponds 
to about 0.75 A. 

It is more illuminating to inquire what is the most probable value of r, 
when 7, is given. For function | of the last section this quantity is 
independent of r,, there being then no radial correlation. It appears that 
for function 2, given by (5.3), the most probable value of 7, is greater or 
less than a certain value 7, according as 7, is less or greater than this 
quantity. For the values of Z and Z’ given above for function 2 the 
value of 7) proves to be about 0.34. The most probable value of r, 
ranges from about 0.34 A to 0.23 A as r, changes from small to large values. 

To examine angular correlation, we must return to the complete 
distribution function (5.1). The importance of angular correlation for 
given values of 7, and r, may be assessed quantitatively in several ways. 
One way would be to compare the probabilities per unit solid angle of 
the two configurations (71, 72, 0) and (r,, 72,7), that is the configurations 
in which the two electrons lie on the same and opposite sides of the 
nucleus along a given radius. The ratio of these two probabilities is 
given by 


2 pPt+2 L bniby 
Pesta) cote S tate (5.8) 
Pitta) 5 p42 3 (1) bad, 
1=0 m>n 


If we neglect all functions 4, for />1, and then use the approximation 
to p+, cos 6,, given by function 3 of paragraph 4, with 


Uiscor ome Mere creeks) (019) 
then (5,8) becomes 
P(ry, Ta) 9) Bs {i ee aren hee (5 10) 
P(r, 12,7) —pr/foJ ” 


where 
/ ae OE (71 +12)/2 


bilbo=* sa Zr, @—4r—4'n 


=4ie(s?@—t2) eo 4-2-2 sech {(Z'—Z)t}, . . (6.11) 
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x being a (negative) constant. For a given value of s, the expression 
(5.11) is always numerically largest when t=0, that is when r)=7. In 
other words, angular correlation is of greatest importance when the two 
electrons are at the same radii, a result which is not surprising. 

If we put t=0, we find that (5.11) has a maximum when 


ry=re=2/(Z"—Z—Z'). . . . . « « (5.12) 

This is the radius at which angular correlation is most important. Using 
the results of Taylor and Parr for the Z’s, this is about 0.674. 

Alternatively, we can estimate the importance of angular correlation 


by calculating the mean value of cos 6,, given the two radial distances 
ry, and7r,. This is given by 


ik COS O;2P(11, 12; 442) 815 
i P(1,, 72; 942) dO» 


Zernerra m/z ore - (48) 


(cos D1), Pian 


8 


If we omit all 4, with />1, we get 


(COS O19). r= 3% 0%1/ (Ho? + 341”). i aaa, (eee Cee 


To a good approximation, the denominator can be replaced by yp? and 
(5.15) then reduces to 3x,/%9, and so is proportional to (5.11). It is of 
interest to calculate the numerical value of this quantity for two special 
cases. In the absence of spatial correlation it would, of course, be zero 
and the average value of 0,, would be 90°. 

(i) If r,=rg=2/(2"—Z—Z’), corresponding to both electrons being at 
the radius (0.67A) at which angular correlation is most important, then 
(COS 6; 9),,, »,= —90.1073 (or 6,.=96° approximately). 

(ui) If r,=r,=2/(Z+Z’) corresponding to both electrons being at the 
most probable radius (about 0.3A), then (cos @,.),, ,,=—0.0687 (or 
6,2=94° approximately). id 

Finally, we may calculate the mean value of cos 6,, averaged over the 
whole of configuration space. This is given by 


cos 0,5= i i, [208 ®12P (rv 815) dr, dry d0,, 


© 39141) 2 ce 
= 2 ore ip I, r2rebabi, ar, dr, . (5.15) 


With the same approximations as before, this is estimated to be 
—0.0665 corresponding to an angle of 93.8°.* 


ee eee eee 


* We are indebted to Mr. Bhattacharya and to Mr. A. Brickstock for some of 
the numerical values given in this section. 


Spatial Correlation of Electrons in Atoms and Molecules 591 


§6. APPLICATION TO ATOMIC AND MoLECULAR STRUCTURE 


In this paper we have been primarily concerned with the spatial 
correlation of two electrons resulting from their electrostatic repulsion. 
It has already been noted in the introduction, however, that for systems 
with more than two electrons, it is also necessary to take into account the 
spatial correlation that arises from the exclusion principle. The effect 
of the latter has already been considered for a number of atoms and 
molecules (Lennard-Jones 1952), 

The next step in the theory must clearly be a comparison of the relative 
importance of the two causes of correlation in a system in which both 
operate. Thus, for the neon atom, the exclusion principle tends, to a 
first approximation, to keep the four outer electrons of the same spin at 
the four corners of a regular tetrahedron (Lennard-Jones 1952). If the 
exclusion principle were the only factor influencing the spatial correlation, 
the two tetrahedra so defined would move independently of one another. 
But if the correlation due to electrostatic repulsion is allowed for, there 
will be a tendency for the two tetrahedra to keep apart, and it seems 
likely that the most probable configuration of the eight outer electrons 
will be that in which each is situated at the corner of a cube. The 
calculations of this paper do indicate, however, that the effect of 
electrostatic repulsion on angular correlation is comparatively small, 
so the relative motion of the two tetrahedra may be much freer than that 
of the electrons within the tetrahedron. It is hoped to publish further 
work on this problem later. 

There are other important applications of these ideas in the field of 
molecular structure. In the paper dealing with the effect of the exclusion 
principle (Lennard-Jones 1952), it was shown that three electrons of the 
same spin in a triple bond, as in the nitrogen and acetylene molecules, 
are most probably in azimuthal planes inclined to each other at angles 
of 120°. If we ignore spatial correlation due to electrostatic repulsion 
(as we do if we use a single determinant wave function), the two systems 
of planes thus defined rotate about the axis of the bond independently. 
Only further development of calculations of the type discussed in this 
paper can lead to further knowledge of the extent to which the six 
electrons are kept in azimuthal planes inclined to each other at angles 
of 60°. q 
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LIV. CORRESPONDENCE 


A New Nuclear-Research Emulsion 
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NucuLeaR-RESEARCH emulsions are often used to investigate the 
properties of the heavy nuclei in the cosmic radiation. Several ways of 
identifying these particles have been proposed and used ; they involve 
a determination of the rate of energy loss combined with that of the 
residual range or scattering. If a particle is known to be moving with 
relativistic velocity (v/c~1), its rate of energy loss is proportional to the 
square of its charge, and identification is then possible by the 
measurement of energy loss alone. 

In most investigations, emulsion sensitive to all ionizing particles is 
needed ; as a result the tracks of heavy nuclei are continuous and it is 
impossible to determine the rate of energy loss directly by grain counting. 
This difficulty has been overcome in two ways. Firstly, it is possible to 
measure the density of delta-rays along the track, a quantity which is 
closely related to the total ionization (Bradt and Peters 1948, Hoang 
1950, Sorensen 1949), or the widening of the core of the track produced 
by low-energy delta-rays (Blau, Rudin and Lindenbaum 1950, Ceccarelli 
and Zorn 1952). The second method is to interleave the sensitive plates 
with very insensitive ones (Bradt and Peters 1950). Grain counting is 
then possible on the tracks of heavy nuclei, and the rate of energy loss 
can be found directly. The advantages of this method are that grain 
counting is simple and quick, that the statistical accuracy for a given 
length of track is greater than that of delta-ray counts, and that grain 
counting is not as subjective a measurement as delta-ray counting so 
that the establishment and maintenance of elaborate conventions is 
unnecessary. The disadvantage of grain counting is that every batch 
of plates has to be calibrated separately. Other disadvantages of the 
method of grain counting as used by Bradt and Peters (1950) are that the 
insensitive plates have to be developed separately by special methods 
in order to obtain the required low sensitivity, and that the tracks have 
to be traced through a stack of plates. 


Fare ee eee eR ee NE as I es ee te 


*Communicated by Sir George Thomson, F.R.S. 


Correspondence 593) 


To overcome the disadvantages of the grain-counting method, an 
experimental emulsion, ‘ Type GO’, has been produced. It differs from 
G5 emulsion in only one characteristic : its sensitivity. 

GO emulsion is processed by exactly the same methods as G5, and it 
is therefore possible to use thick layers, and to make sandwiches of 
GO and G5 emulsion on the same glass base. As can be seen from fig. 1 
(Plate X XVII*), the interface between the two emulsions is clearly defined. 

The relation between grain density and energy loss in a batch of 
GO emulsion has been determined by grain counting in a sandwich 
plate consisting of 4004 G5 emulsion covered with a 200-u layer of 
GO emulsion. It is shown in fig. 2. The plate was developed with 
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Grain density as function of energy loss in GO emulsion. The errors are 
standard deviations assuming a Poisson distribution. The grain densities 
to be expected for various relativistic nuclei are indicated. 


amidol-bisulphite (Herz 1952) to an average minimum grain density 
(in G5) of 26 grains/100 . The tracks used in the calibration were those 
of alpha-particles and singly-charged particles ending in the emulsion, 
taking all singly-charged particles to be protons. 

GO emulsion has been designed in order to simplify the work of 
identifying heavy nuclei, to increase the accuracy for a given track 
length, and to increase the emulsion-to-glass ratio in the stacks of plates. 
It appears to satisfy all these requirements. Moreover, in addition to 


*Hor Plate see end of issue. 
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its application in cosmic-ray research, it may well be useful in nuclear 
physics, because it is an emulsion of low sensitivity which can be used 
in very thick layers, either alone, or sandwiched with layers of electron- 
sensitive emulsion. Work is continuing on materials of this type and it 
may be found desirable to produce emulsions of somewhat different 
characteristics. 
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LV. Notices of New Books and Periodicals received 


Linear Computation. By Paut 8. Dwyur. (John Wylie & Sons Inc.) 


Tsis book deals with the solution of sets of linear algebraic equations and with 
related problems such as the calculation of the inverse of a matrix and the 
formulation of the characteristic equation. It is concerned almost exclusively 
with ‘ direct ’ methods, i.e., methods equivalent to eliminating the variables in 
turn, and iterative or other indirect methods are referred to only briefly. The 
book starts with a discussion of the approximate nature of the numbers used in 
most numerical calculations and the author defines five or six kinds of approxi- 
mate number. While this discussion may be useful in clarifying ideas it cannot 
really be said to lead anywhere since once the author gets down to business he 
follows the general practice of using rounded-off numbers in which the error is 
unspecified and usually unknown. The book gives a minute account of the 
different methods with their variants which have been proposed for the solution 
of sets of linear equations. This involves much detail since differences in, for 
example, the order in which the multiplications or divisions are performed can 
make all the difference between a good and a bad numerical procedure. 

The author’s own interest in the subject began in the field of statistics and he 
explains in the preface that the book is designed for readers with very little 
mathematical equipment. Determinants and matrices are not used in the 
earlier parts of the book but they are defined and their properties are developed 
ab initio in the latter part. The reviewer cannot help feeling that those who 
wish to go into the subject in the detail offered in this book could be expected to 
have more mathematical background and would have preferred a more concise 
treatment. Such a treatment would also have been more useful to people 
needing practical advice about how to solve a set of equations ; at present such 
people are likely to be bewildered by the amount of material in the book, and 
they will have to search to find the author’s conclusions about the relative 
advantages of the various methods he gives. The book is well provided with 
references and exercises for the reader to work. M. V. W. 


Dielectric Breakdown of Solids. By S. Wuirenzapv. [Pp. 271.] (Oxford 
University Press.) Price 25s. 


Tus is a most valuable book on dielectric breakdown. It contains a great 
amount of material of considerable interest to the engineer and the experimental 
physicist working in this field. There are six chapters respectively on method of 
treatment, intrinsic breakdown, thermal breakdown, breakdown caused by 
discharges, electro-chemical deterioration and on dielectric breakdown in 
practice. The reader interested in the fundamental side of this field will welcome 
the inclusion in Chapter II of a certain amount of theory, mainly based on 
Froéhlich’s work. Theory has accounted for a number of features of breakdown 
in a qualitative way and in some cases reasonable quantitative agreement with 
experimental data had been found. Nevertheless further, well controlled 
experiments, in particular on time effects of pre-stressing on intrinsic breakdown 
are wanted in order to obtain a more complete picture. The great difficulties 
in the necessary experimental technique as well as in the interpretation of the 
results will certainly be appreciated by readers of this monograph. 

Finally, it should be mentioned that the material in the present volume on 
dielectric breakdown is almost completely different from the contents of the 
earlier (1932) book on the same subject by the same author. soe 
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X-Ray Analysis of Crystals. By J. M. Bisvonr, N. H. KoLKMEWER and C. H. 
MacaitLavry. (Butterworths Scientific Publications.) Price 50s. 


Tus volume is the first English edition of a well known Dutch work published in 
1938 and revised in 1948. For descriptive convenience only, its subject material 
may be most easily compared with ‘ The Crystalline State—Vol. I’ by W. L. 
Bragg but differs chiefly in respect of its simpler and less detailed approach, and 
by way of its greater emphasis on structural crystallography ; also the more 
recent material and references which it contains. As such it is a useful compre- 
hensive text-book for the science student whose course includes crystallography, 
and also for the research beginner. 

Starting with a simple introduction to x-ray diffraction directions, and 
crystal symmetry, it leads to the intensity of diffraction structure determination 
by trial and error and Fourier methods. In order to keep the text as lucid as 
possible, some of the details of the subject, such as the concept of reciprocal 
space, line- broadening etc., are presented in a fairly lengthy appendix. 

The second, slightly shorter portion of the book, is perhaps the most successful 
and is devoted to an excellently illustrated explanation of basic structure types. 
After a brief classification by bond-type, and a discussion of simple structure and 
the concept of atomic size, the dihalides, silicates and long chain compounds are 
emphasized as examples. 

It is difficult to keep a comprehensive account of any subject up to date, and 
it is encouraging to meet reference to modern researches such as dislocation in 
metals, haemoglobin and neutron diffraction to quote just a few. 

Although the text is reasonably complete in itself, the more serious worker is 
presented with twenty or so references at the end of each chapter, a fair pro- 
portion of which are later than 1945. This ample documentation makes the 
book interesting to the more serious student of the subject. 


[The Editors do not hold themselves responsible for the views 
expressed by their correspondents. ] 
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Growth spirals on the (111) faces of gold crystals obtained by precipitation. 
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Appearance of Tracks in G5/GO Sandwich Emulsion. @ and 6 are heavy nuclei, 
c is a radio-active star. The scale of 6 and c is the same, and different 


from that of a. 
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